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SYNOPSIS 

_ i'piPPUSIOn MOAMYTIO EEAGTION 111 A POES" A Thesis 
Suhmixted in Partial Piilfilment of the Eequirement s for the 
Degree of Doctor of Philosophy in Ohemical Engineering hy 
ji.B, Lulkarni to the De-par tnent of Ghenical Engineering, Indian 
Institute of leciinology, Kanpur, August 1972. 

A study on ”Dif fusion and Konoatalytic Eeaction in a 
Pore" has heen undertaken with a view to smalyse a chemical 
reaction in a porous mass in terms of that in a single pore. 

The shrinking pore approach illustrates the role of physical 
factors in heterogeneous kinetics with special reference to 
the possible inverse effect of temperature on the overall rate 
of a nonsolid forming reaction. The scope and limitations of 
single pore analysis in relation to the existing methods of 
investigation are discussed in the introductory note. 

The revievf of literature enumerates some of the important 
developments in the field of heterogeneous fluid -rsolid reactions 
v/ith special regard to micromechanistic, phenomenological and 
empirical aspects. The subject matter of graphite oxidation 
is reviewed with the objective of arriving at a suitable 
esqoerimental procedure and appropriate mathematical spproxi^^ 
mations. Some of the controversies regarding the mechanism 
and order of carbon gasification, role of mass transport, 
catalytic activity of impurities and control regimes are 
highlighted. Important aerospace and nuclear applications of 
graphite are also included in the review. 
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Expo iiinient al efforts are oonfined essentially to weight 
loss measurements during the course of oxidation of graphite 
hy air or oxygen, h plane surface of graphite is allowed to 
react in a stagnant air column and the resulting wei^t loss 
curves are used to determine the apparent first order rate cons- 
tant. Eor studying the reaction in a pore, the external surface 
of a hollow graphite cylinder coated vfith a protective layer 
and the internal surface is allowed to react. Difficulties 
encountered therein are discussed. Thermo-G-ravimetric ilnalyser 
(T.G.h,), Veeco vacuum evaporator and gas chromatograph are used 
in the experimentation. 

One dimensional diffusion equation with the reacting flux 
boundary condition is solved to obtain a suitable mathematical 
algorithem for the determination of surface rqte constant k, ' 

Phis k is used in building various mathematical models to des- 
cribe the reaction in a pore. After exhausting the possibility 
of analytical and seminumerical solutions, implicit and explicit 
finite difference approximations for the relevant set of sim-ul- 
taneous partial differential equations are tested,' Von Ueimaann ^s 
stability criterion for explicit schemes and 1^ error estimates 
for Milton Lees*, Rose’s and Modified Ba.okward difference schemes 
are discussed in detail. The manor portion of the thesis work 
is concerned with the numerical solution of a n-umber of stable 
finite difference equations and with the minimization of machine 
time on IBM -7044 and IBM-560, Conclusions are then drawn by 
critical analysis of the numerical data on pore wall profile 
and axial concentration distribution. 







1, INTEODUCTIOU 


Even though the parameters pressure and flow invariably 
have a ■telling effect on the rate of a he'terogeneous fluid-solid 
reaction, temperature seems to have occupied a more prominent 
position in the scientific and technological investigations of 
a number of industrially important reactions such as reduction 
of iron ores, roasting and smelting of metal sulphides, com- 
bustion of solid fuels and solid propellant, gasification of 
coal or oil shale, regeneration of catalysts, oxidation of 
iron by steam, calcino.tion of carbonates, flourination of 
uranium ores, chlorination of titanium dioxide, and aluminium 
oxide, sulphation of carbon and carbonization of calcium oxide* 
Besides the inherent technical difficulties associated with 
high pressure and high flow systems, the factors responsible ' 
for this undue importance of ■ temperature seem to be a lack of 
qualitative appreciation and quantitative understanding of 
the role of mass transport in a porous reacting solid , system. 
‘Traditionally the relative -importance of physical and chemical 
factors have been studied on one or more of the following lines; 

(a) Langmuir type of adsorption isotherms, and 
their application to mechanistic and kinetic 
s-fcudies , 

(h) Extension of well established method of catalytic 
kinetics, • 

(c) Microstudies of reactions involving single crystals. 
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(d) Phenomenological approach to shrinking core 
models, 

man the rates of diffusion through fluid film and 
porous solid are both very fastj the overall rate of a fluid- 
solid reaction is solely controlled by the inherent chemical 
reactivity of the solid reactant. Based on the langmuir type 
of adsorption isotherm, it is possible quantitatively to 
describe the mechanism of uncatalysed heterogeneous fluid- 
solid reactions under this condition. However, depending on 
the mechanism, the resultant rate equation involves more than 
two arbitrary constants, some times as many as seven-. In 
selecting these constants for each mechanism, the curve repre- 
senting the rate equation of the -favoured mechanism which best 
fits the experimental data is chosen. Because of the unavoid - 
able intrinsic scatter of the experimental data, it is clear 
that in many instances little meaning can be attributed to the 
magnitude of the adsorption equilibrium constants, the frequency 
factor, and the apparent activation energy. Often, the difference 
in fit may be so small that it is very difficult to determine 
v/hether it is simply due to experimental error or truly due 
to the difference in mechanism. Purthermore, an a,ltemative 
mechanism may fit the data equally well necessitating an - 
additional extensive experimental sea,rch for the correct 
mechanism. Although a correct mechanism which is difficult 
to determine, will allow extrapolation to .conditions not 
actually investigated, there is no reason why simple rate 



eq.uations which, fit the data satisfactorily should not be 
used provided no extrapolation is done beyond the range 
investigated. Based on the adsorption isotherms, the order 
of fluid-solid reactions can be shown to vary from zero to 
two depending on whether the gas reactant is strongly, modera- 
tely or weakly adsorbed. Experimental studies also indicate 
many reactions to have a similar range of the order of 
reaction depending on conditions such as structure of the 
solid, reactivity, pressure, temperature etc.[1] ,[2 ] , The 
very fact that adsorption is only one of the six steps involved 
in the rate process of a heterogeneous reaction, any attempt 
to describe the entire process kinetics by means of just one 
step is bound to have serious inadequacies. Moreover, the 
assumptions made in the Langmuir approach are too ideal to be 
applicable to an engineering fluid-solid system. Thus adsorption 
isotherms can at best bo useful in mechanistic studies. 

Although comprehensive theories of catalytic reactions 
are available, studies on noncatalytic reactions have not been 
advanced and organized into a distinct and coherent branch of 
chemical reaction engineering. This slow development stems 
from littlo understanding of intrusion of complex physical 
phenomena on chemical reactivity on one hand and the hetero- • 
geneous surface properties and pore size distributions on the 
other hand. One of the major difficulties involved in the 
analysis of noncatalytic systems is the tine varying moving 
boundary condition and the resulting "geometric instability" [3] 
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of the porous structure. Some attempts have been made in 
recent years to develop the concept of effectiveness factor 
and Thiele parameter for simple first order irreversible 
fluid-solid reactions ?;ith no change in number of moles during 
the course of reaction [3 ,4 ,5 ,6]. However, too many 

simplifying assumptions involved in such analyses restrict the 
application of these developments to a number of important 
designs. Unlike catalytic reactions, most of the reactions on 
the solid surface cannot be truly reversible, since the pheno- 
menon of solid deposition in the reverse reaction does not 
necessarily provide the solid product with a structure that 
is exactly the same as that of the original solid reactant, 

YOiile for the first order reactions it is possible to provide 
on additional term to account for the ’’equilibrium hindrance”[ 31 
the same may become very complicated for other orders of reaction-^ 
The thermodynamic and transport properties of a noncatalytic 
system are dependent not only on temperature, pressure, and 
composition, but also on the time varying pore structure. The 
limited applicability of Wheeler's single pore model [7] or 
Smith's macro-micro pore .diffusion model [ 8] is also due to the 
timevarying pore structure. Thus it is felt that the study 
of variation of pore structure in general and a single pore 
in particular, during the course of a fluid -solid reaction, ! 

may become a key factor in understanding a heterogeneous reactive 
system. It is in this background the present work has been 


undertaken. 
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In the recent past gas-solid reactions haTe been studied 
on microscale by means of a number of new methods such as 
etchpit expansion [93 , measurement of varying thermo e.m.f. 

w ■ 

at the reacting junction [l0] , micro cinographic technique [nlj 
molecular beams, olectronmicroscopy, EER, ESR techniques etc. 
While these studies yield results directly applicable in the 
field of solid state chemistry, the.ir utility in obtaining 
design data, for macro-heterogeneous systems appears quite 
doubtful and remote. Until well developed statistical theories 
to predict the macrodata from microdata are made available, 
the results of the above micro investigations can be considered 
only of potential future use. 

A good number of engineering investigations have made 
use of a "shrinking core model" to obtain the kinetic data 
required in the design of . heterogeneous reactors [ 12 , 15 , 14 ■, 
15, 16, 17, 3], Eespite a wide application of the shrinking 
core model, it is by no means the panacea for complex non- 
catalytic reactions. One of the most common assumptions made 
in this model is the "quasi-steady state of shrinkage" which 
may not be valid for systems with high pressure and low solid 
reactant concentration. When the solid contains enough voidage 
or develops enough voidage during the course of reaction, so 
t;|Lat the diffusion into the interior of the solid is free and 
f^st, the shrinking core model will no longer be applicable. 

This model cannot be used uniformly when shifting of the rate 
controlling regimes takes place. In such cases Wen [53 tias 
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considered the reaction to occur homogeneously throughout 
the solid. He has developed a generalized model which also 
has its limitations to the extent that the rate controlling 
regimes are replaced by "stages of reaction" i.e, using 
different sets of equations and solutions for initial, inter- 
mediate and final stages of reaction. The shrinking core model 
assumes a sharp demarcation between the reacted and unreacted 
zones and invariably uses a heterogeneous rate constant defined 
on the basis of surfa.ce reaction only. However, when "geome- 
trical instability" appears, the reaction interface is no longer 
smooth, and the mathematical model may at best become an approxi- 
mation of the real process. If the results obtained for a 
regular shaped pellet are to be generalised and applied for 
an irregular shaped pellet, the manner in which the timevarying 
pore geometry changes v/ith the core geometry (spherical, 
cylindrical or rectangtar) poses a problem in the shrinking 
core model. Thus it appears that some of the limitations of the 
shrinking core model may be overcome by a quantitative under- 
standing of the variation of pore surface area, porosity,, and 
the poresize distribution during the course -of fluid-solid 
reaction#, 

* 

After having realized the importance of a "shrinking 

pore model" in making the "shrinking core model" more versatile 

and useful, and after realizing the fact that a single pore 

forms the most natural fundamental unit of an entire solid core, 

one may think in terms of two lines of approach for the 
— : — — — ^ — ^ — ^ — — 

^The term hereafter relhhe 'shrinkage of lengthy 



7 


development of a ”s3arinking pore model", 

(a) llie deterministic approach which starts with the 
study of a single pore and then tries to describe the variation 
of the entire pore structure in terms of its fundamental uniti 

(b) The stochastic approach which makes use of time- 
dependent probability distribution parameters to describe the 
variation of a pore structure which consists of a large number 
of pores of various sizes and shapes with a complex interlinkage 
among them. 

In the present work, an investigation on a single pore 
is undertaken with the objective of making a modest beginning 
in the deterministic approach so that a basis for a new method 
is made available for future research investigations. The 
essence of the present effort is to determine both the axial 
concentration and the radius of the inner surface of the pore 
as functions of time and axial distance when a fluid-solid 
reaction is taking place within a pore. These functional 
relations in turn a.re used to predict some of the design data 
pertaining to heterogeneous kinetics. Attempts are made to 
substantiate the numerical calculations with the experimental 
investigations. Some of the areas in which the above single 
pore analysis may be used are indicated below. 

The Kozeny Carman equation for the pressure drop across 
a packed bed or a porous medium can be derived either from a 
"grain model" or from a "pore model". It is surprising that 
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similar exaiistive treatment has not been extended to reactive 
porons media. Only recently Szekley and Evans [is] have 

r 

reported such developments for the reduction of ironioxide 
by carbon monoxide. It is interesting to note that the grain 
and pore models are implicitly synonymous with the volumetric and 
the surface reaction rate constants respectively. As long as 
one of the products of reaction is a solid, (i.e. when the 
porosity during the course of reaction remains more or less 
constant) it is not of great consGq.uence whether the voidage is 
represented in terms of interspace betw'een the grains or the 
volume of pores. However, when a fluid-solid reaction produces 
only gaseous matter, pore structure undergoes a continuous 
change, and pore model - at least intuitively - may have an 
edge over the grain model, because -physically it is a pore 
which forms the fundamental unit of a reacting porous solid* 

It should hence be logical to follow the changes occurring in 
a single pore and' then extend this knowledge to a multiple pore 
system by making suitable estimates of the mean pore dimensions, 
their distributions and their linkage. Alternatively, single 
pore data may enable us to predict the rate constant in a porous 
medium from the known rate constant in a nonporous medium. 

Except the new experimental technique of Westenberg 
and Walker [ 19 ] j there probably exists no accurate method for 
the determination of gaseous diffusivities at temperatures 
above 600°K, Frequently, experiments to determine diffusivities 
are conducted in nonreacting gaseous systems with constant 
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compositions, and the data is used for the reacting systems. 

If the first order surface rate constant is known, the method 
developed in the present investigation can he used to determine 
the gaseous diffusivities in a multicomponent reacting mixture. 
This method acquires a special significance v^hen it is desired 
to consider the dependence of diffusivities on concentration, 
and on the possible dissociation of gases at high temperatures. 
Secondly, if the gaseous diffusivity through a porous reactive 
solid medium is to be determined, one needs to know the change 
in cross-sectional area available for the diffusion as the 
reaction progresses. The model presented in this work may serve 
as a basis for the estimation of changes in cross-sectional 
area of a reactive porous solid at various stages of reaction. 

With the advent of high speed computers and newer finite 
difference schemes for the quasilinear parabolic differential 
equa.tlon, the time-consuming trial and error procedure of 
numerical calculations should not be too problematic,. At the 
same time, the experimental procedure is simple and can be adapted 
easily to new systems. 

As a corollary, one can utilize the same analysis to 
determine the first order surface reaction rate constant, if 
the diffusivity and rest of the parameters are' a.lready known. 

The conventional experimental methods for determining the 
heterogeneous rate constants employ a system wherein the 
reacting fluid either flows over the solid or remains stagnant 
on it. If it is a flow system, the true rate constant is 



inYaxiably an extrapolated value at Infinite velocities* If 
the gas is stagnant, it is rather difficult to eliniinate the 
convection currents. To avoid both these difficulties, the 
present method employs the inner surface of a pore of a smaller 
diameter vjhcra the radial symmetry for temperatLire and concen- 
tration more or less frees the system from the effects of 
convection; and hence., the hulk diffusion terms do not appear 
in the mathematical model* 

Whenever diffusion and chemical reaction occur simulta- 
neously,- the transition zone hotv/een the diffusion controlled 
rogime and the chemical reaction controlled regime becomes an 
important factor in the design of a reactor. If the experiments 
are carried out at various temperatures, the present method may 
serve as a starting point in distinguishing the diffusion 
controlled regime from the reaction controlled regime. As will 
be seen in the conclusions., steeper profiles of the reacting 
pore-wall would indicate broadly a diffusion controlled regime 
find the flat profiles, a reaction rate controlled regime. It 
will also be seen that the conventional method of deciding the 
regimes by the criterion of an almost constant reaction rate 
equalling the diffusion flux at the pore mouth, may at times 
give misleading conclusion. The present method may serve as 
an additiona,! check to resolve the controversies if any regarding 
the nature of the transition zone. 

By introducing some extra terms for natural convection, 
or bulk mass transport with or without flow, the present model 



can "be modified fo predict the rate of chemical corrosion of 
tube v/alls, storage tanks, or the metallurgical laddies or 
any vessel which contains a corrosive fluid. Of course, the 
mathematical complexity of the problem increases v\?hen we try 
to incorporate the terms for convection or flow, Nevertheless, 
the modol with its numerical solution provides a fundamental 
system theory approach to corrosion problem. This may also 
help to isolate the effect of corrosion from that of errosion 
in metallurgical laddies or runways. So far, the methods for 
predicting the corrosion of reacting walls are mostly empirical 
and with tho present model, a beginning in the fundamental 
approach can be made. 

Sometimes, the specific reaction rate constants depend 
on the geometry and the heterogeneity of the sample surface [20, 
2l 3 , This nay be due to the difference in surface potentials 
f;)r differont geometries, or may be due to the difference in the 
crystallite orientation on the machined surface of the solids, 
Quito often, corners act as the points of stress concentration 
and hence are more susceptible to chemical reactions'. If such 
a phenomenon is anticipated or observed, one can cast or 
machine a hollow cylinder out of the solid and carry out the 
reaction in the pore so that variation of rate constant with 
geometrical configuration can be detected or confirmed, for 
example, the orientation of carbon atoms on the surface of a 
machined graphite depends on the relative geometry of the sample 
Therefore, the surface rate constants in the initial stages of 
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oxidation of rectangular, spherical, and cylindrical graphite 
saniplos are observea to be slightly different froni one 
another [ 2] . 

One of the important considerations for deciding the 
optimum particle size is the ratio of the overall rate of 
reaction on the inner surface of the pores to that on the outer 
surface of the solid, A mathematical model with its analytical 
solution presented in this work deals with this ratio which 
also indicates the ratio of mass transfer due to convection on 
the exterior surface to that due to molecoilar diffusion inside 
the pore. Variations in porosity during the course of reaction 
can also be studied in the present experimental set up. 

In view of all the possible interesting applications 
of the study of noncatalytic reactions in pores, one may think 
in terms of a voolid-liquid system or a solid-gas system for. 
the experimental investigations. The latter is preferred in 
this v/ork bocausc of the following reasons: 

The magnitude of gaseous diffusivities is generally 
larger than that of liquid diffusivities and this helps to 
ostahlish appreciable axial concentration gradients during the 
course of reaction. Besides, the larger magnitude of gaseous 
diffusivities provide a good spread for temperature as a 
variable in the numerica,! calculations. Pressure can be a 
meaningful variable in case of a gas-solid system while this 
possibility is almost absent in a liquid-solid system. Heats 
of wetting for liquids in general are larger than those for 
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gasos3 and lionco the changes resulting exclusively from the 
\ 'irlcit iono ill iihe heats of wetting tend to become relatively 
less sigiiiu-f icant, -The formation of a curvature on the liquid 
surir.ee inside a pore and the effect of this curvature on the ■ 
equilihriu.m vapor pressures are very typical of a solid-liquid 
system while these are almost absent in case of noncondensable 
gases in Cf ntact with solids. In fact almost all physical 
phenomena resulting from heterogeneity of the surface of the 
solid become loss significant in case of a gaseous system. In 
a restricted sense (i.e, confining to a few systems which are 
under present consideration) it may be noted that the nonideal 
behaviour with respect to activity coefficients is niore severe 
in case of liquids than in case of gases, for gases one 
can have better control on the para,metors like flov;, temperature., 
p.ressurG and composition, low viscosity of ga.ses reduces the 
momentum entry length at the^pore mouth so that the end effects 
due to flow are minimized. All these factors tend to reduce 
the side effects in the main stream of ana.lysis and experimen- 
tation so that the certainty and the validity of the conclusions 
are enhanced. After having decided upon a gas-solid system, 
the particula.r reaction of oxidation of graphite by air has 
been selected for a number of reasons. 

Well defined crystal structure with minimum number of 
defects is an important criterion for the choice of a solid, 

Next to diamond, graphite is the most ideal carboniferous \ 
material for the study of reactions of carbon with oxygen. The 
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structure of gropliite is ■well defined and is ea.sily controlled 
at its foruativo stages, and therefore one ma^ be confident 
in building the atonic or molecular models for studying the 
surface properties. Extensive studies pertaining to the kinetics 
and nochanisn of graphite oxidation are already available and 
therefore one can be more certain about the validity of the 
experimentally detcrciincd parameters like order of i:eaction^ 
primary products of oxidation, specific reaction rate constant, • 
and activation energy. This, in turn, will naturally minimize 
the uncertainty about the main stream of analysis,. Also prior 
knowledge concerning the transition temperature zone between 
the diffusion controlled regime, and reaction rate controlled 
regime, pr.”)vides some additional guide lines for simplifying 
the oxpori mental and numerical procedures,- Substantial knov^ledge 
about the variation of properties of graphite with respect to 
por-.r.ity is alrocdy available [22] and this supports some of the 

r 

assuriptions made in the mathematical modelling. Graphite has 
a good machinable property and hence samples of the desired 
shape and size can easily be obtained. Air happens to be a 
versatile n;iterial for which the thermodynamic and transport 
properties are tabulated over a wide range of temperature and 
pressure. This again enhances the certainty of the main stream 
of analysis. Besides, experimental and analytical simplicity, 
the subject of oxidation of porous graphite in recent days is 
gaining considerable importance in the field of nuclear power 
generation and rocket propulsion where graphite is used as a 
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material of construction. And lastly the analytical techniques 
iisoO ill f ollowin.i''': the course of chemical kinetics of graphite 
oxidation are faix’ly simple irrespectiye of "whether one uses 
weight laoo measurement technique and/or gas composition 
measurcmani; teclmiquo. However the aforementioned considerations 
should not in any way undermine the contraversy that exists in 
the field of graphite oxidation [ 2,2^’. 

n * 

Before concluding this introductory. note, it is felt 
necessary to mention some of the limitations of the approach 
introduced in this vrark. The most importan-t aspect in this 
respect is the specific meaning of the' title "diffusion and 
noncatalytic reaction in a pore". The term "diffusion" in 
general v/ould moan mass transfer under concentration gradient 
v/hich may ho resultant of a partial pressure gradient, total 
pressure gradient, temperature gradient, velocity gradient, 
Durfa.ee potential gradient, ionic, molecular and steeric 
factor gradients. However in the context of the present work, 
mass transfer d\ie to partial pressure gradient alone will he 
considered as diffusion. The diffusion controlled regime 
mentioned in this work_ refers to this aspect of diffusion only. 
The reaction is carried out at a constant pressure and it is 
assumed that the number of moles do not change during the course 
of reaction and therefore the total pressure gradient cannot 
possibly exist in the system. Mass transport v/ith unequal 
absolute pressures at the ends of the transport path, often 
called permeability, is specifically excluded from consideration. 
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In fact, the nse of a macropore for the experiment should 
rjle out the above possibility* It goes ■without saying that 
fciio Knudsen diffusivity is also out of consideration* It may 
DO noted that this microdiffusivity is proportioncS,! to the 
a.bsoluto temperature to the poY/er of half, while the ordinary 
diffusivity is proportional to a power variously reported as 
from about 1,5 to 2 , 0 , Tlie assumption of an isothermal 
condition of reaction, eliminates the possibility of thermal 
diffusion. In absence of external fields, and negligible ■ 
dissociation of molecular oxygen below 1.500°C, other types of 
diffusivities are not likely to be significant* According to 
IJichols [ 24 ] inspite of the above limitations, we may still 
be cnoountorlng two types of diffusivities, The first one is 
ordinary, bulk, or the Poiseulle diffusion and the second one 
is surfOiCa diffusion, a process of mobile adsorption which can 
occur in parallel w'ith cither of the above types of diffusion. 
It is independent of total pressure but is proportional to a 
complicated function of temperature, Wicke and Kallenbach [25} 
claim that their data on diffusion of carbondioxide through 
charcoal at atmospheric pressure shows a sizeable fraction of 
such- surf ace diffusion in the total mass transport. Surface 
diffusion is generally considered to occur mainly at a 
temperature close to the boiling point of the diffusing gas or 
in extremely fine pores at other temperatures [ 24 ], In this 
work, neither of these conditions is prevalent, and therefore 
surfa.ee diffusion is less likely to interfere with the kinetic 
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data. Before concludinguit may be worth recalling another 
kind, of diffusion wnich occuisv;hen molecular dimensions are 
'.ipproximatoly oquqI to the dia,nieter of the pore# This kind 
of activated diffusion is exponentially dependent on temperature. 
This also is a_ remote possibility on the open surface of the 
graphite walls. It may be noted that only gas phase axial 
diffusivity is considered in the present model, For reasons 
mentioned in the chapter on mathematical modelling, no radial 
diffusion in the gas or so'lid phase is considered. 

The term "chemical reaction" refers to a simple first 
order isothermal irreversible surface reaction without any change 
in the number of moles, lo complex or compound reactions due 
to secondary changes or inhibitive poisoning actions are dealt 
with. The appfirent density of the solid during the course of 
roaction is assumed constant even though occasionally there is 
,a denoity variation in the direction of "depth" from the surface 
of the solid. The "equilibrium hindrance factor"i33 due to 

r 

the change in the pore structure or in the geometrical stability 
of the solid does not enter the picture hecause of the assumption 
of irreversibility. The term "pore" in the context of this 
analysis refers to a macro pore only, A perfect cylindrical 
geometry with two open mouths is taken into account at the 
start of reaction. 

Even through reference is made to a number of potential 
applications of the present study, the actual work is restricted 
only to the prediction of wall profiles and concentxatioh 
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pi’oftlcs in a {:raphite cylinder, the inner surface of ?i?hich 
is ccn tiniicnsly reacting with air. The author is well aware 
of tne linitations experienced in representing the surface 
reaction hy a single first order rate constant at all chrono- 
logical stages of reaction - more so at higher temperatures. 

Of c '^urse, the inherent inadequacies associated with the con- 
cept of the order of a heterogeneous reaction and Pick's 
second law of unsteady state diffusion are common to most of 
the invuc-tigritior.n jxi lieterogcneous systems, 

Tims the first chapter of the thesis introduces the 
single p<oro method with its scope and limitations. Chapter 
t'wo gives a brief review of literature on heterogeneous gas-,, 
solid systems in general and the specific reaction of graphite 
oxidation by air in particular. The third chapter describes 
the oxporiraental methods both for the determination of the ., ■ 

firot order rate constant and the study of the reaction in a 
single p'-re. Chapter four presents various possible mathematical 
models and the assumptions involved there in. The fifth and 
the final chfipter deals with the resrilts of computation from 
which some conclusions are drawn, .. : 


* * * 



2. CI' LITERATURE 


Depending on the extent of stress that is laid on 
tn.e tLeinnodynamic , mechanistic, kinetic, heat or mass transfer 
pr?.ncipLes, the existing literature on heterogeneous fluid- 
solid reactions may be classified into three broad categories s 

(a) IviGchanistic studies dealing essentially v/ith the micro 
aspects. 


(b) Empirical and semiempirical studies with utilitarian 
objectives, 

(c) riienomeno logic hi stii.diGs taking a maciro approach to the 
process of slmult.aneous diffusion and chemical reaction. 

Uhilc microstudios are quite useful in understanding the 
mochsnisci!, of reo.ction, heterogeneity of the solid surface and 
tho inadequacies associated with the micro -experimental 
teoiniqu'-o a.ro two factors which limit the scope of the micro- 
rioch'inistic investigations. Even if accurate microdata is 
available, its application to macrosystems poses a number of 
pliysical and riathonatical problems for which no satisfactory 
solutions exist. The details of the microtechniques and the 


difficulties involved in their utilization for design purposes 
are discussod in the later part of the review. 


As regards tho empirical and semiemplrical data .,, they 
are too numerous to be comprehended. Every investigation has 
a limited objective of its o?m and invariably scientific 
generalizations are not possible, Eor exa.mple, thin films of 
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lud rmd ^vrarhito act as bearing lubricants and their oxidation 
1- air or oxygon [26, 2?] has some unique features which are 
not ccnmon to heterogeneous fluid-solid reactions in general. 

The data, of Bowden and Young [26] regarding the effect of 
degasification of grapnito surface on its lubricating properties 
cannot bo generalised and used elsev/here, Maberry et.al, [ 28] 
iuivu Lstuiiou tho porf or nance of gr 3 .phite nozzles used in rockets 
atui their data on the corrosion and errosion of internal wall 
suriace ax’o valid only oyer a restricted range of operating 


conditions. Mast of ohe inyostiga,tions pertaining to the 
reduction %;i iron '■.re pellets and the sintering effects involved 
therein are by and laipo empirical in nature. The data regard- 
ing the change in 'xorous structure of the solids [ 29 ] during 
tlio CDuruc of heterogeneous reactions is still %t large empirical< 
The qiuuo titative effect of irradiation on the oxidation chara- 
ctoriJtics af a. graphite moderator [ 30] has usually been 
described by merxjs of empirical formulae. In spite of enormous 
efforts, the study of the effect of impurities and crystal 
defects on the course of a -heterogeneous reaction has remained 
an empirical scienco. Many more instances of empiricism can 
be cited from the litera-ture, but compilation of such data is 
not the a,in of this reviev;, . 

*'ConceT)tual modelling'* or a ’’phenomenological formulation” 
has been an essential theme of a large number of engineering 
and technological investigations. In fact this has been the 
most versatile approach to the synthesis of design data. The 
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rocoGc.ra'.v nr- rhonr: t.ioal "lodel is built by talcing into consi- 
der"^,! "11 I* i-lo ^ phanoKcnr!. or steps tba,t are involved in 
t.'U,' alien ;.i' r.-oc + ioit. Owing to the close resemblance between 
catalytio .'uvl noncntalytic systems, the sclicDie suggested by 
Hougnn ond Watson [51] for a catalytic reaction has served as 
t:io basis for fornulating a. parallel scheme for noncatalytic 
rcactiiins. Won [3] has described a surface reaction to consist 
of tl;,e f . llnwing stops: 

(l) DilTus'i jH of the fluid reactants across the fluid-film 


surrouridlr.g tho solid, (2) Dif 


id, (2) Diffusion of the fluid reactants 


tlircu'rh p'ur- us solid layer of reactant and/or product, (5) 


a,dscrption of the fluid rca.ctants at the solid reactant surface, 
(4) ol'emical reaction v;ith tho solid surface, (5) desorption of 
tlio fluid products from the solid reaction surface^ and (6) 
diffusion of toe product away from the reaction surface through 


the portiua c 'lid modi a and through the fluid film surrounding 
tho solid. Since fcho above steps take place consecutively if 


any one of thorn is much slower than all others, that step 
becomos rate— do tor mining. In other words, for a given set of 
reaction conditions, that step which offers the highest: resistance 
contrcDls the overall rate. In phenomenological formulations, 
it has been a custjm<ary practice to combine the steps 3 and 5 
v/ith the step 4 so that the surface reaction can be represented 
by means of a single equation involving only iwo constants, y 
Thus tho reactions may be external mass transfer controlled, 
inter or intra particle diffusion controlled, chemical reaction 
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controlled, lieat transfer controlled or mixed controlled. 

However single-step rate determining processes are limiting 
cases. The majority of fluid-solid reactions are influenced 
simultaneoiisly ty moro tiian one step, j 

Thus a rigorous treatment seems unattainable even for 
the solid of the simplest geometry, under isothermal conditions 
of roaction. Partly this is due to the intricate relationship 
among the rates of chemical reactions and the rates of mass 
and energy transfer. Overall rates are influenced not only - 
by the rate of chemical reactions occurring in or at the surface 
of a solid and by the mass transfer rates of fluids through 
the solid as well as across the fluid-film surrounding the 
solid, but also by factors such as solid reactivity of crystallite 
oriento-tion, crystallite size, surface characteristics, impuri- 
ties etc. Besides a great number of difficult problems exist 
in pro,ctical systems, such as the changing size and shape of 
the solid during the reaction and formationoFa product around 
the solid reactant which may crack or ablate or change the 
porestructure drastically. In addition, the complex velocity 
profile of the surrounding fluid makes the problems of mass and 
energy transfers to the solid reactant more difficult to analyse. 
Often one encounters the problems of adulterating the chemical 
reaction constants v/ith the structural parameters and other 
physical factors. If during the course of rea.ctton there is 
an appreciable change in the solid density or if there is a 
structural nommiformity or a thermal stress in the solid, the 
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solid particles may undergo cracking, ablation, and peeling 
off, during tlie course of reaction. Most of the reversible 
fluid-solid reactions cannot be truely considered as reversible, 
since the phenomena of solid deposition in the reverse reaction 
do not necessarily provide the solid product with a structure 
that is exactly the same as that of the initial solid reactant 
in the forward reaction. This is seen in a number of fluid- . 
solid reactions which have been known t'o exceed their normal 
thermodynamic equilibrium [ 32 , 33 , 34 ]. If this is the state of 
affairs for a simple rea,ction system, one can well imagine the 
magnitude and complexity of the situation when compound and 
complex reactions with nonisotherraal and nonisobaric conditions 
are provailent in the system. Therefore it should not be ' 
surprising if there has been a lack of co-ordination from the 
phenomenological approach resulting in, for the most part, frag- 
mentary end some times obscure nature of theories of non-catalytic 
reactions vjhich have not been advanced and organized into a ‘ 
distinct and coherent branch of chemical reaction engineering* 
Although there 'are a number of diverse and complex 
fluid-solid reaction systems, ?/en[3] has shown that it is 
possible to group them phenomenologically in two ways. One 
grouping is based on the phase combination of the reactants 
and the products, and the other is based on the manner by which 
a reaction progresses. 
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Glassification Tsased on the phases of reactants and products ; 
(a) Solid reactants fluid products 

Some of tile examples for type A are pyrolysis of carloonaceous 
materials, combustion of double-base propellants, and thermal 
decomposition of some organic or inorganic compound especially 
explosives such as ammonium nitrate, 

^ + SHgOCg) 

(B) Solid reactants »■ fluid and solid products 

Examples are 

CaCO^(s) ^ CaO(s) + COgCg) 

CaSO^-lHgOCs) ^ CaSO^*^ H20 (s) +|H20(g), 

CaOClgCs) — -»■ CaCl2(s) + ^ 02(g). 

(C) Fluid and solid reactants *■ fluid products 

Examples are 

C(s) + S2(g) ^CS2(g), 

+ C(s) — ^ — -^NaClT(l) +H2(g), ' 

2Al205(s) + 6Cl2(g) — 4AlCl^(g) + 302(g). 

.(D) Fluid and solid reactants : solid products 

e,g, GaC2(s) 4- ■— *" »- CaCir2(s) + C(s) 

(E) Fluid and solid reactants fluid and solid products.: 

Examples sire 

Fe20^(3) + 300(g) 2Fe(s) + 3C02(g) 

2 ZnS(s) + 302(g) — ^ 2 ZnO(s) + 2 S 02 (g) 

4EeS2(s) rf 11 02 (g) 2Fe20'^(s) + 8 SOgCg) 
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Classification "based on the control rep:imes; 

(A) Diffusion controlled reactions or heterogeneous surface 
reactions: ?/lien the chemical reaction rate is very rapid and 
relatively the diffusion is sufficiently slow, the zone of 
rea.ction or the chemical reaction propogating front is narrowly 
confined to the interface between the unreacted solid reactant 
and the product. Also when the porosity of the unreacted solid 
is very small so that the solid is practically impervious to 
the fluid reactants, the reactions will occur at the surface 

of the solid or at the interface between unreacted solid and 
the porous product layer. Such a reaction may be , considered a 
heterogeneous surface reaction or the diffusion controlled 
reaction. Since there exists a sharp demarcation between thd 
reacted and unreacted zones, an unreacted shrinking core model 
is conveniently used to describe the reaction* Yagi and Kunni's 
ano.lysis [35] is a classical exa^mple of a shrinking core model, 

(B) Reactions controlled by the chemical reaction rates or the 
homogeneous reactions: In many cases, the solid contains 
enough voidage to permit free passage of the fluid reactant and 
the fluid product and the solid reactants are distributed homo- 
geneously through out the phase. In other words, the effective 
diffusivity of each fluid component through the entire solid 
phase is sufficiently large and also the rate of mass transfer 
to the periferi of the solid is sufficiently fast* Ihen under 
such circumstances it may be rcaBonable to assume that the 
reactions between fluid and solid components are occurring 



26 


homogeneously tlirougla out the solid phase. Homogeneous model 
which defines the rate constant on the basis of an apparent 
volume of the solid, may he employed for the mathematical 
analysis of such systems, Levenspiel’s model [ 36] has been a 
significant development in .this direction, 

(C) Reactions controlled by heat transfer rate or those 
accompanied by phase changes of solid components or by evolution 
of volatiles; In some cases, because of a rapid exothermic 
reaction and poor heat conduction, a phase change of the solid 
component, tales place prior to the chemical reaction, The solid- 
reactant simply melts or sublimes before it is brought into 
contact with fluid reactants, Ihe phase change may proceed with 
pyrolysis or devolatilization which is often associated with 
gasification and combustion of solid fossil fuels. In such . 
cases a homogeneous vapour phase reaction takes place either 
around the periferi of the solid or in the solid product layer 
formed around the unreacted solid core, Harsimhan’s heat 
transfer model [37] which fitted the experimental data of 
Satterfield and Reakes [38], is a good example of such systems 
.controlled exclusively by heat transfer rate. 

(D) Reactions governed by mixed controls or intermediary 
controls; In a strict sense, no heterogeneous reaction belongs 
exclusively to the types A, B or 0, This is because no real 
system can possibly have diffusivity, thermal conductivity or 
rate constant of zero or infinite value. Therefore, the 
categories A, B, and C can only be mathematical approximations 
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to the real situation. On a molecular scale, the distrihution 
of solid reactant in the solid phase cannot be considered 
homogeneous. At best, it can be considered as an ensemble of 
sma.ll lumps of reactant distributed throughout the solid phase. 
The overall reaction rate depends on a number of parameters 
such o,s the intrinsic reaction rate constants and orders, the 
sorption characteristics and the transport properties of the 
fluid reactajots in and around the porous solid^ the structure 
of the solid, and the distribution of the small lumps of solid 
reactant. As has already been pointed out, a master model 
which incorporates all these factors seems unattainable. In 
fact: most of the cases where the experimental data do not 
fit into any one of the above three models belong to the 
category h. To gain bettor insight into such situations, the 
recent trend has been not only to isolate kinetic parameters 
from diffusion parameters, but also various types of diffusion 
parameters from one another, ©bkarn and Doraiswamy [59] have 
cited a number of examples of mixed controls in an illustrious 
way, 

Insplte of a. number of serious short comings of the 
above classification, almost all phenomenological investigations 
have relied on it. Each reaction is studied separately at low, 
intermediate, and high temperatures so that separate analyses 
for chemical reaction rate contio lied regime j mixed control 
regime, and diffusion controlled regime can be made use of. 

The presence or absence of an "ash layer" makes a distinct 
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phenomenological difference and thus fhe classification based 
on the phases of reactants and products becomes implicit in 
the mathematical model, 

Gokam and Doraiswamy [39] have systematicallj presented 
the significant developments in the study of ash forming hetero- 
geneous fluid-solid reactions. In most of these reactions, the 
reacted solid retains its original geometry and dimensions, and 
the reaction front recedes from the outer surface to the centre. 
Such topochemical reactions under the quasisteady state conditions 
are characterized by a linear rate of movement of the reaction 
zone, A critical reviev of some shortcomings of the shrinking 
core model used in those reactions, has been presented in a , 
number of recent publications [3, 18, 40, 41]. Existence of 
a sharp interface between the reacted and unreacted regions is 
a fundamental assumption of the shrinking core model. In 
practice, the demarcation is never sharp. It is not uncommon 
to find a hemispherical core shrinking into an irregular shaped 
cone during the course of reaction. A circular ring demarcating 
the reacted and unreacted zones during the course of reaction 
may become a wavy ring of any odd shape. Thin line of demar- 
cation may grow into a diffuse structure. With such developm.ents, 
surface reaction rate constant will have little meaning. Wen [ 3] 
has given a quantitative estimate of the accuracy of a pseudo- 
steady-state approximation which may not be valid for a high 
pressure and low solid reactant concentration system.^ In other 
words, the initial unsteady state time need not alwa^-^s be 
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negligible compared to the pseudo-steady-state time of reaction. 
In the case of an unreacted-core shrinking particle, when the 
rate per unit area decreases as the conversion becomes greater^ 
the nonuniformity of structure will tend to be smoothed out 
because the rate of reaction is less at a point of deeper pene- 
tration, On the other hand, when the rate of reaction per unit 
area increases as the penetration is deepened, greater uneven- 
ness of the reaction surface will result. This is termed the 
"geometrical instability" and ?/as first pointed otit by Cannon 
and Denbigh [42], Shrinking core model does not normally account 
for such instability. The assumption of Gq.uimolal coimtercurrent 
diffusion through the ash layer or the constancy of effective 
binary diffusivity may not be valid if there is a change in 
the number of moles dLiring the course of reaction or if the 
concentration of diffusing species is high. At times the bulk 
flow due to pressure gradient cannot be overlooked in the 
shrinking core model. Moreover the irregular surface boundary 
of the exterior solid core sets in complex velocity profiles 
in the adjoining fluid layer. This will result in unpredictable 
mass transfer coefficients at the fluid-solid interface* 

Shrinking core model must be used with caution when shifting of 
rate controlling regimes takes pla,ce, For example, a fluid- 
solid reaction originally under the diffusion-controlled -regime, 
may shift to the chemical-roaction-controlled regime resulting 
in an overall rate, many folds lower than the -original. In 
such cases the unreacted-core-shrinking model, of course 
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represents only a part of the reaction history. The reacting 
porous solid is considered to he semi- infinite with a macro- 
scopically unidirectional movement of the reaction zone i.e, 
one dimensional diffusivity. But in practice the propagation 
of reaction in a porous solid may be. multidirectional. Iso- 
thermal conditions, uniform pore structure etc, are assumed 
after careful consideration of the individual system. No genera- 
lizations are possible in this respect. 

The concept of the order of heterogeneous reactions. -— 
though empirical - has been a powerful tool in describing the 
most complex phenomena in terms of relatively simple expressions. 
It has been aptly pointed out by Arts and Mehta [ 45] that "t^«gh-‘; 
the concept of reaction order is at best a difficult one to 
justify - and this is particularly true in the context of 
heterogeneous catalysis - its undoubted success in organizing 
data on complex reaction systems has led to its wide spread 
use”. Thus zeroth order reaction, in which the rate is indepen- 
dent of concentration, may be approached in a heterogeneous 
reaction when the surface is almost completely coveredi while at 
low coverage the reaction may be of the first order, fractional 
orders have been found to fit experimental data with tolerable 
accuracy. Negative orders have been used to represent inhibitory 
effects. At times an "apparent" order of reaction is used to 
represent all sorts of physical and chemical factors and 
possibly "lumping" in complex reaction schemes. The surface 
rate constant and activation energy which are closely associated 
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witli the concept of order of reaction, invariably are adnl- 

^ * 

terated with pore diffusioni formally, the order of- a surface 
reaction is based on the concentration of reacting fluid in 
the gas phase. Hovrever, the concentration of solids too has 
been incorporated in the rate expression by Wen [3]. This will 
have, some advantages in the construction of a homogeneous model. 

At times structural parameters like solid surface heterogeneity, 
imeven poresize distribution, varying porosity and crystallite 
orientation are included in the empirical rate eq.uation but 
satisfactory solutions to such problems s,re yet to be obtained. 
Whenever a heterogeneous reaction tahes place in a mixed 
controlled regime, the concentration gradient of the reacting' • 
fluid along the depth of the solid is significant. In such 
situations, attempts have been made to define and develop the . 
concept of effectiveness factor for none atalytic. reactions, 
hue to some intrinsic differences in the solid reactivity and- 
the diffusion characteristics, the effectiveness of the effective- | 
ness factor for a noncatalytic reaction has not been as signi- ; 
ficant as that for a catalytic reaction. Mathematical complexity | 
increases as the changes in pore structure during the course of 
reaction become more severe. And the most severe changes perhaps -i 
are experienced in a reaction with no solid product. Depending . ; 
on the conditions of a reaction, the basis for the definition 
of the effectiveness factor may be the void volume of the, entire 
solid or the void volume of the unreacted mass only, .Similarly 
the maximum value of 'the fluid concentration used in the definition 
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of the effectiveness factor may vary from system to system* 

One may consider the concentration at the interface betv^een 
the reac.ted and unreactod mass, or the concentration at the 
fluid-solid : interface, or the hulk concentration in the fluid 
phase. Time dependant pore surface area adds one more dimension 
to the complexity of the problem, le?/is and Paynter [44] have 
shown that for most of the catalytic batch reactors, the 
difference between the steady state and the unsteady state 
effectiveness factors is less than ten percent. The same is 
not true for all none ataly tic reactions. Therefore one has to 
identify the common bases and conditions before making any 
comparative evaluation of reactive systems in terms of their 
effectiveness factors,. 

Por the first time, the classical work of Thiele [45] was 
extended to noncatalytic reactions (gG,sification of carbon by 
c.arbondioxide) by Peto3:-sen [46], ?/alker [ 2] gave a formal 
definition of tho "Thiele Utilization factor" for the first 
order irreversible reactions of carbon and tabulated the 
expressions for sphericeil, cylindrical, and plane blocks in 
various temperature zones, Weisz and Prater [4?] developed a 
general criterion for the importance of pore diffusion with 
chemical reaction in porous solids, 4nd this is one of the 
most useful practical results of the concept of effectiveness 
factor for first order reactions, Petersen [48] has pointed 
out that for other order reactions, particularly for adsorption 
mechanisms, the criterion can be drastically different. 
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Petersen [49] fuctlicr gave an asjmpto-tic solution for the 
generalized effectiveness factor and a generalized criterion 
for the importfmce of pore diffusion, Aris[4] then shov;ed 
hovv' this asymptotic Thiele parameter can be normalised so that 
the effectiveness factor is asymptotically proportional to the 
inverse of the Thiele parameter when its value is large., 

Bischoff [50] has arrived at a still more general criterion of 
Weisz Prater type but not v/ithout the presence of an unevaluated 
integral in the criterion equation. All such investigations 
deal with general situations which are common to both catalytic 
and noncatalytic systems. Recently Ishid a and Wen [6] have 
developed the concept of nonisothermal effectiveness factors 
exclusively for the noncatalytic reactions, Prom the effective- 
ness factor plots, three types of instabilities are identified; 
the geometrical instability, the thermal instability?’ due to 
metala,ble' temperature., and the instability resulting from 
abrupt shiftings of the rate controlling regimes, Calvelo and 
Cummingham [ ^ have presented a similar generalization of the ' 
effectiveness factor for fluid-solid reactions representing the 
moving bo’undary, the continuous layer and the moving layer 
conditions. However, suitable modifications of such generali- 
zations are yet to be made available for thq particular reaction 
of graphite oxidation. 

Existing literati.ire on graphite oxidation runs invariably 
in two distinctly divisible strea^ms, the first one dealing 
essentially wi’th thero'idynamics, kinetics and mechanism of 
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oxidation, and the second one confining itself to the process 
of diffusion of oxidant through porous graphite, OnLj recently 
some attempts have been made to combine these two streams to 
illustrate the effect of simultaneous diffusion and chemical 
reciction on the continuously changing internal wall profile 
of a single pore in particular, and on the entire pore structure 
of graphite in general, [ 29, 50, 51, 52], Accordingly,, the 
current review will also follow a similar line i,e. to begin 
with the tv/o stroarris will be dealt with separately, and then an 
attempt will be made to combine the two, 

The controversies and difficulties associated with the 
process of oxidation of carbon and carbonaceous materials are 
dealt with exhaustively by a number of investigators [ 53, 54, 55, 
56, 57, 58, 59] . In spite of these difficulties and may be 
because of them, the investigations on gas carbon reactions 
are really extensive and voluminous extending back to the middle 
of the nineteenth century[7]. It is neither possible nor 
necessary to include all the extensive work in a short review* 

So stress will be la,id only on these investigations that are 
relevant to the present work. Rather, an effort will be 
made, to consider only notable advances in this field'. Ihe 
most extensive review 'on the subject of Chemistry and Physics 
of carbon has been made available by Wa Ike r[ 60a, 60b, 60c, 60d* 
60e], In the following part of the review, the thermodynamics 
of the relevant carbon gasification reactions prececle the 
kinetics and mechanism of such reactions. 
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gHBMOMMICS 

Oxidation reactions in general are exothermic and are 
accompanied by highly negative changes in free energy. Graphite 
oxidation is no exception to this and as such the thermodynamic 
eq.'uilibrl'um is not lilcely to be a crucial factor in evaluating ■ 
the progress of reaction. However, thornodynamic calculations 
are necessary to choose proper ranges of operating pressures 
and temperatures so that the product of oxidation is essentially 
carbondioxide. By doing so, it is possible to neglect the rest 
of the side reactions, and assume a relatively simple kinetic 
model for further computations. . 


TABLE 2-1 


Ho, 

Reaction 

h'l 

i K.Cal. 

AG K.Cal 

1* 

C(s) + 02 ( 8 ):^ 

11^002 (g) 

■r94.P3 

-94.260 

2, 

C(s) +»0o(g).y— 

i:rco(g) 

-26,62 ■ 

r32.808 

3. 

C0(g)+i02(g)“ 

irC02(g) 

-67.41 

-61,452 ' 

4. 

c(s) + co 2 (g);:i: 

=:;r2co(g) 

+40.79 

+28.644 

5. 

C(s) + H20(g)“ 

ZZ!-C0(g)+H2(g) 

+31.14 

1:2 1,820 

6, 

C0(g)+ H20(g>?zz 

i::2-C02(g)+H2(g) 

-09,65 

-06.817 


Eor a system consisting of carbon, oxygen and water 
vapour, various possible primary and secondary chemical reactions 
with their respective heats and free energy changes at 18°C and 
,1 atmosphere are listed in Table Ho. 2-1, In stochio metric 
chemical equations, carbon is assumed to be in 3 form with 
zero heat of formation. Based on this, various forms of amor- 
phous carbon are reported to have positive heats of formation 
ranging from 1,7 to 2, ,6 K.Oal/mole [45] which can be neglected 
in the preliminary calculations. 
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Independence of reactions can be worked out by writing 
tlio above set of reactions in the following forms . 


c 

+ 

Op 

- 

CO 5 


z 

4 

z 

4 

z 

= % 

c 

+ 

P2 

4* 

z 2 

— 

CO 

4 

z 

4 

z 

= Z 

z 

+ 

tOo 


CO 5 

4 

CO 

4 

z 

4 

rr 

£j 

= z 

c 

+ 

ir 

4 

CO 5 

— 

2C0 

4 

z 


z 

= z 

z 

+ 

z 

- 

oo| 

4 

GO 

4 

HgO 


F 

"‘2 

= z 


h.B. Z indicates a zero in the above equations, 

Ihcrofore the matrix of stoichiometric coefficients is 

1 1' -1 "O 0 0 

1 i 0 -1 0 0 

0 i -1 -1 0 0 

10 1-211 

.1 0 0 -1 1 -1 

0 0 -1 + 11-1 

By going through the standard procedLire suggested by Aris [61] 

we obtain the m 0 ,trlxs 

11-10 0 0 
0 1 -2 2 0 0 

0 0 0 0 0 0 

0 0 -0 0 0 -0 

0 0-1 1 1-1 

0 0-1 1 '1 -1 

Obviously the number of independent reactions is TWO, Here 
the reactions 1, and 5, may be considered to occur independently, 
Brom Table 2-1 it is clear that at 18°C and 1 atm., the 
reactions 4 and 5 are thermodynamically not feasible while the 
reactions 1,2,3 and 6 are feasible, Bor some of the reactions free 
energy changes are nearly zero, and hence it becomes necessary to 
investigate the effect of temperature and pressure on the 
equilibria of these reactions.. Variation of the equilibrium 
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constant witli respect to temperature is indicated in 


Figure 2-1, From this one con easily see that the reactions 



of graphite with carbon 
dioxide and with water 
are not likely to be 
significant within the 
range of temperature 
under consideration. 
Further, reaction 1 has 
an equilibrium constant 
which is much higher 
than that of any other 
reaction, Therefore, 
from a limited thermo- 
dynamic point of view, 
complete combustion to 
carbon dioxide appears 
to be more favourable. 


From the data of reference [62] it is clear that the oxidation 
of c^arbon to carbon monoxide and ^arbon dioxide is not restricted 
signif icsntly by equilibrium considerations even at temperatures 
upto 4000‘^K, It may also be noted that temperatures above 
2000®K should be, avoided if reactions with carbon monoxide 
and/or water are to be suppressed. 


Again, from table 2-1 ■ it can easily be seen that the 
effect of pressure on the equilibrium constant can be studied 




Equilibrium CO/CO 2 ratio as a Equilibrium Product Steam Ratios 

function of temperature and for Reactions 5,6,4. 

pressure for the reaction 

C + OO 2 = 2G0 


EIGI3RE 2-.2 


EIGIIRE 2-3 
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only in the case of reactions 2, 3, 4 and 5 where the number 
of moles does not remain constant during the course of 
reaction. Over the range of temperature unvder consideration, 
the fugacity coefficienis of the gases, up to a pressure of 
100 atmospheres, remain very close to unity and hence the 
effect of pressure is not likely to be significant, From 
actual calculations also we find that the numerical values of 
K do not vary significantly v;ith pressures even upto 1000 

Jr * 

atmospheres. However, the variations in equilibrium compositions 
of carbon monoxide and carbon dioxide are indicated in Figures 
2-2 and 2-3 [2]-. From Fig. 2-3 it is seen that at temperatures 
above 1200°K and at atmospheric pressure the conversion of 
carbon dioxide to carbon monoxide by reaction 4 is unrestricted 
by equilibrium considerations.' At elevated pressures the 
possible conversion markedly decreases; hence, high pressure is 
likely to suppress this reaction and promote reaction 1. For 
the carbon steam reaction, it is seen in Fig, 2-3 that the ■■ 
amounts of carbon monoxide and hydrogen, which can be produced 
above 1100°E and upto a pressure of 100 atmospheres, are 
essentially equal, even when the possible, side reactions are 
considered. However, as in the carbon-carbon dioxide reaction, 
the possible extent of conversion of steam to carbon monoxide 
and hydrogen decreases with the increasing total pressure. 

Since the moisture content of the system used in the present 

investigation is very low, and the steam-carbon reaction is 

* 

less likely to occur under noncatalytic conditions, the presence 


4-0 


of hydrogen in the gaseous mixture is ruled out. Consequently, 
the reaction' of methane formation is not included in the 
thermodynamic review. 

After considering the effect of temperature and , pressure, 
one may be tempted to go in for the analysis of a complex 
equilibrium for the reactions 1,2,3 and 4 to arrive at the 
equilibrium compositions of carbonmonoxide, carbondioxide and 
osygen in the multicomponent reaction mixture. But uncertainty 
of kinetic constants and the inherent coupled diffusion effects i’* 
associated with porous graphite render the complex analysis less 
meaningful. This may perhaps be the reason why the analysis 
of complex equilibr'lum is not to be found in the vast published 
literature on this topic. 

Since the sorption of gases on the surface of graphite 
is the primary step - a.tleast from a chemist’s point of view - 
in all graphite oxidation reactions, the thermodynamics of 
chemical reaction remains incomplete vhthout the thermodynamics 
of sorption. The interaction between the graphite surface and 
the oxygen molecule can be any one of the following; 

(a) monolayer adsorption 

(b) multilayer adsorption or condensation 

(c) chemisorption in presence or absence of other 
impurities 

(d) chemical reaction producing volatile products 
Energetics of each of these interactions will naturally be 
different from one another. 
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While tliere exists an alimdant literature [ 22, 26, 63> 

64 ] on adsorption of gases on carton tlacks, Yery little is 
reported about the adsorption of gases on graphite [65] • This 
should be expected because incomplete graphitization has been 
proved to produce active polar groups on the surfaces of 
active carbon particles [23]. Correlation of adsorption pover 
and catalytic activity of various finely divided carbons shows 
that these lie in the sequence; active charcoal > soot > coke 
or retort graphite, for near ideal graphites, the heat of 
adsorption has been estimated as fev; ergs per gramme mole of 
oxygen, and no calorimetric experiments have been possible in 
this direction. So, f'^r all practical purposes, the adsorption 
type of interaction on gra.phite surface can safely be neglected,. 

Coming to the next stronger type of interaction of oxygen 
with graphite, much controversy is found, not only about the 
physicochemical nature of interaction, but also about the 
nomenclature used to describe this interaction [64]-. The 
current review ' folio v/s Ubbelohde's nomenclature [ 23 ] to 
describe these interactions which can be classified into two 
broad categories, . 

(i) Graphite Oxide which is the product of chemisorption 
of oxygen on graphite in presence of hydrogen, which may be 
represented by a general reaction like xC + yO + zH'“— ^ CxOyHz, 
(ii) Surface oxide which is the product of chemisorption 
of oxygen on graphite in the absence of hydrogen, which may 
be represented by a general reaction like xC + yO— 0^0^, 



WithoTJit going into the controversies regarding the 

precise chemical nature of graphite oxide or surface oxide, 

it may he worthwhile to recollect some specific instances . 

where the heats of interaction of graphite and oxygen are 

indicated. A graphon which provides a homogeneous hydrophobic 

2 ~ 

surface of 100 metre /gram is supposed to have a medium surface 

energy of 110 ergs/cm [60b page 221] . Heats of immersion in 

,ca.se of HaDBS solution are of the order of 2,5 cal/gm* of the 

same graphon [60b page 219] , With H-butanol, in. aqueous 

" " * 2 

solutions, it is of the order of 20-80 ergs/cm. Bull H.J,. and 
M.H* Hall[66] in a very nudimentary way calculated from bond . ■ 

energy considerations the heat of chemisorption of oxygen as 
112 K*cal/grammo mole, which was in reasonable agreement with 
the experimental values, R, Nelson Smith [67,] has reported 
the variation of the heat of chemisorption with temperature 
and the amount of oxygen chemisorbed. The first quantities of 
oxygen admitted to carbon surface at room temperature are chemi- 
sorbed with exceedingly high heats of formation (70-100 Kcal/ 
mole) of the same order of - magnitude as the heat of formation of 
carbondioxide. These heats decrease rapidly as additional 
oxygen is added, finally levelling out to about 55-60 Kcal/molej 
when physical adsorption sets in.. The heat of formation of ^ 
complexes formed at 400°C (calculated from heats of combustion) 
is 32.5 Kcal/gram atom of oxygeni independent of the amount 
of oxygen combined. All these figures indicate the obvious 
diversity of thermodynamic data on sorption. Engineering 
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studios have therefore- preferred to group together the 
steps of sorption and desorption with chemical reaction, 
thus enabling the usage of a single equation for the entire 
phenomena, 

lE IGIiMISM MD KlhBTICS 

It is needless to say that the doubtful nature of the 
chemical composition of the surface complex, surface hetero- 
geneity and surface impurities affects the mechanism as much 
as it does affect the thermodynamics of sorption. However,, 
a broad general scheme of a chain reaction mechanism can be 
put forv/ard in terms of a "first step" and a series of 
"succession steps" [68, 69, 70 , 7 l] , Many interesting suggestions 
regarding "edge attack", "surface attack", foreign element 
attack at a "claw" or a "hole", breaking off of volatile 
molecules such as carbon monoxide, chemisorption of atomic 
oxygen or carbon monoxide, and formation of lamellar type of 
compounds, have preceded the development of modem concepts 
of the nature of defects in solids and the different ways in 
which foreign atoms can be incorporated in solids. In view 
of this and because of lack of precision in characterising 
the defect state of carbon or graphite, the situation regarding 
the "first step" and the "succession steps" can be considered 
open for further investigations. 

Upto i960, the maior question conceming the mechanism of 
carbon-oxygen reaction has been whether carbon dioxide is 
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a primary product of reaction or a secondary product resulting 
from the gaspiiase oxidation of carbon monoxide. By now, 
however, a large majority of workers have agreed that both 
carbon monoxide and ca-rbon dioxide are primary products, 
and the primary carbon monoxide-car bo ndioxide ratio increases 
with increasing reaction temperature. If high velocities of 
oxygen are used to avoid stagnation due to molecular diffusion, 
partial pressure of oxygen does not significantly affect the 
primary carbon monoxide-carbondioxide ratio, A number of 
correlations are aVcCilable for the variation of the primary 
carbon monoxide - carbon dioxide ratio with temperature, but none 
of them is really versatile because it has not yet been esta- 
blished that the magnitude of this ratio is solely a function 
of temperature and is independent of -the carbon reacted and 
the impurities present in it. It is interesting to note that 
the evidence for the formation of a carbon suboxide ( 0 ^ 02 ) 
as a first volatile product at present refers to graphite with 
considerable radiation damage, and may apply only to highly 
defective carbon net works [54, 72] « Another suboxide (C 2 O) 
has also been proposed in oxidation mechanisms particularly 
at low temperatures [73], The mechanism of carbon-carbon dioxide 
reaction, which is itself not understood properly, cannot 
possibly be used to understand the nature of the primary 
product of oxidation. 

Before proceeding with the details of kinetics, it 
may be useful to list some experimental techniques used in 
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the study of graphite oxidation. One of the commonLy employed 
methods of following the kinetics of reaction is to observe 
pressure changes as a function of time [74], Since the 
products are volatile, another convenient way of measuring 
reaction rate is from the loss of weight of the graphite 
sample [ 64, 75, 55] , as recorded by moans of a sensitive 
microbalance operable at variable temoerature, pressure and 
gas composition. Recently developed- instruments provide a 
gas analyser to continuously record the composition of out- 
going gases, hot only can this give information regarding the 
carbonmonoxide-carbon dioxide ratio, but it also provides a 
double check on v/eight loss measurements. Very high velocities 
of gas have been utilized to eliminate the effects of diffusion 
while sta,gnant gas is \3sed to introduce the effect of molecular 
diffusion only* Inhibitors like phosphorous trichloride are 
employed in the gas stream to arrest secondary reactions* 
Radioactive tracers have also been used to follow the 
reaction [7^] . 

A number of classical teebniques have been employed to 
isolate the effect of catalyst from the crystal defect. An 
experimental method applied in high temperature studies is to 
use a carbon filament as a specimen which can be heated by 
the passage of an electric current and the temperature can he 
measured by means of an optical pyrometer [ 77] » This not only 
solves the problem of material of construction at high 
temperatures, but also eliminates the parameter of porosity 
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firoffi experimentation wlien very thin filaments are employed,- 
EectmiTly [60a] attempts have heen made to record, the , .c hang es 
in t?aormo-e.m.f , developed at the gunction of graphite crystal, 
to follow the kinetics of chemisorption of oxygen* Optical 
and electron microscopy techniques have teen very useful in 
studying the etc]ipit expansions during the course of oxidation. 
J.M. ihomas [ GOa] has reviewed important experimental techniques 
in this direction, and has indicated how etchpit expansion 
studies have helped to differentia.te between the rate of oxi- 
dation in planar and perpendicular directions at the molecular '■ 
level.. Recently, [ 1l] , a novel technique of a motion picture 
study of catalysed and uncatalysed reaction between carbon and 
oxygen has been effectively employed i:.i Pennsylvania State 
University, Coming to macro stvidies, huge channels and nozzles 
of graphite are subjected to oxidation studies under flow and 
nonflow conditions to find out- the suitability of graphite as a 
material of construction in nuclear reactors and rockets[52, 28] 
With the available experimental data and the postulated 

mechanism, various attempts have been made to interpret acti- 

& 

vat ton energies, rate constants order of reaction, at different 
temperatures. In absence of diffusion effects and secondary 
reactions, the rate of reaction (weight loss of carbon) can 
be considered the same as the rate of surface rearrangement 
of the carbon-oxygen complex to a ra.pidly desorbable picoduot. 

In principle., in an equimolal gas solid reaction,, the order of 
reaction with respect to the oxygen pressure should vary 



■be'bween ONE and ZERO as ©, tlie fraction' of the surface 
covered hy carhon oxygen complex changes from ZERO to ONE, 

At a constant temperature, 9 varies from ZERO to ONE as the 
pressure is increased from zero to higher values, thus 
causing a corresponding variation of order from ONE to ZERO. 
Similarly the reaction temperature can also affect the 
order. The reaction which is of zero order at low temperatures 
and at a given pressure can become of first order at the . 
same pressure and sufficiently higher temperatures. 

The majority of results [74, IQ, 79, 80, 81, 82, 83-, 84, 
under varied expericiental conditions show the carbon- 
oxygen roaction to be first order or close to first order, 

Erom the previous discussion, it implies that under all the 
experimental conditions used by the above authors, the fraction 
of the total active carbon surface occupied by the oxygen 
surface complex at any given instant during the reaction 
a.p pro aches zero. Almost all the authors have given experimental 
evidence to this effect. Two notable exceptions are found to 
the carbon-oxygen reaction being first order, G-ulbransen and 
Andrew [5^ , working with spectroscopic graphite,, have found 
that, at reaction temperatures of 450* C and 500*0 and at^ 
pressures below 0, t5 cm. of Hg, the order is nearly^ zero. 

They do state further that at pressures above 10 cm, of Hg, 
the reaction is first .order, Blyholder and Eyring [ssl reacted 
extremely thin.' coatings of graphite-^ supported on a ceramic 
base, with oxygen at 800*0 and a pressure less than 100 |j. Hg, 
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Prota liffiited data tiiej have concluded that the reaction 
is of zero order. Andrew’s statement that the reaction 
order increases v/ith increasing pressure is difficult to 
explain on theoretical grounds because, in any heterogeneous 
reaction, .with increasing pressure © tends to OTE and the 
order of reaction to ZERO. However, Blyholder and Eyring 
prescribe a half order reaction for Gulbransen and Andrew’s 
data at 450°C over the entire range of pressure. In the 
present work, low temperatures have been avoided so that a 
first order reaction may be miiformly assumed. 

The oxidation rates for carbon filaments increase with 
tempers,ture upto about 750°-1000®C, Betv/een this range and 
about 1700°C rates remain of first order, but decrease with the 
rise in tempera.turo, Eigure 2-4 illustrates results obtained 
by different authors who have studied the oxidation of carbon 
filaments. As the temperature increases, the rate goes 
through a, maxima? at still higher temperatures the rate increast 
again, but the order now becomes ZERO, Different steps appear 
to exert their control in different ranges of temperature, 

Erom an oxidation rate equation of the form Rate = A + BPo 2 > 
where A and B are constants at a pa,rticular temperature and E 02 
is the partial pressure of oxygen-, a rate controlling acti- 
vation energy of 37.5 Kcal/g.mole has been inferred [53]. 
Alternative interpretations can be that this energy refers 
to the desorption of chemisorbed carbon monoxide molecules, 
without which the succession steps cannot occur, or to the 
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"breaking of carbon bonds as carbon monoxide is separated 

* 

froia the carbon network* A rate expression of the form 
Rate = A y^Po 2 bas’ been proposed as an alternative to fit 
the same experimental data and also other results [85]^ 

This gives an activation energy of about 42 Kcal/mole, This 
form of appa,rent rate expression could originate from the 
diffusion of oxygen into the cracks and pores of the solid, 
as a rate controlling process [ 73 . Attempts to avoid compli- 
cations of porosity by using na.tural graphite of high density 
(2,26 g/cc) still fail to evade the problems of total pore 
surface area and particle size [ 86] , Where diffusion of 
oxygen into graphite controls the rate of oxidation, it has 
been shown that the order of reaction is (n+l)/2 ?/here n is the 
true order of the surface reaction [85] , Under these conditions 
the observed activation energy should be half the activation of 
the actual surface process. When very thin layers of graphite 
are oxidized so as to minimize the effect of diffusion into 
the interior, zero order reaction rates with respect to 
pressure were obtained wdth the calculated activation energy 
of 80 Kcal/mole, Since this value is also found in the zero 
order oxidation rate that prevails above 1700°C for graphite 
filaments, this appears to refer to a fundamental step in 
oxidation such as the breaking of C-C bonds in the surface, 
figure 2-5 illustrates a scheme of activation energies 
for different steps involved in graphite oxidation when the 
rate is controlled solely by resistance to chemical reactivity 
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ejid not ty diffusion. Apparently, step 2 appears tlie slowest, ■ 
However, Rossberg[87] suggests that the slowest step in the 
reaction is the separation of an oxygen atom from the reacting 
species, E© justifies this by comparing the activation 
energy of the reaction C.+ l/2 Og — CO (58 Kcal/g.mole) with 
the dissociation energy associated with the breaking of half 
mole of oxygen into an atom of oxygen i.e. l/2 O 2 — ^ 0 (59 

Kcal/atm of oxygen). This line of argument appears inconsistent 
because it presumes that irrespective of whether & is ZERO or 
ORE, dissociation of oxygen is tie slowest step. As discussed 
before, the implication of the zero order reaction is that the 
overall gar-:if ication rate is being controlled by the rate of 
removal or rearrangement to a desirable product of the surface 
oxygen complex and not by the rate of its formation. Even 
when the reaction is of first order, it is doubtful whether 
Eossberg's concept has any significance, since it would appear 
unsound to compare half the dissociation energy of oxygen 
with the activation energy of reaction. Moreover, the 
dissociation of oxygen just does not figure in any of the 
steps mentioned in the reaction mechanism. Perhaps a more 
reasonable explanation of relative activation energies for 
the reactions between carbon and oxyge-n containing gases, is 
indirectly hinted by Long and Sykes [88], Iheir contehtioh'--' 
that the exothermic ity of interaction of a gas molecule with a - 
carbon free site in case of carbon-oxygen reaction is nearly 
twice that for carbon-steam or carbon-carbon! dibxide'..reac:tions 
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fits v,fell in the pattern of actiyation energies shown in 
higiiro 2-5, 

for some carbons it is contemplated that the reaction 
with carbondioxide proceeds at a measurable rate above 600°C, 
Ho?rover, with eleotrograde 'graphite of density 1.6 gm/cc, 
the reaction with carbondioxide' is "unlikely because excess of 
carbon monoxide molecules in the gas phase compete with 
carbondioxide for active sites, thus poisoning the reaction. 

This is supported by tracer experiments with radiocarbon [ 89] 

But with synthetic graphite diffusion into internal pores could 
lead to a different explanation of the inhibitory effect [ 75] , 

In other words, inhibition is imperative and hence the isolation 
of reaction of oxygen from that of carbon dioxide is needless . 
Reaction with water vapour - as has been pointed out in the 
thermodynamic revie?/ - is a more remote possibility than that 
wfith carbon dioxide, and hence need not be considered. It is 
to be noted that while the presence or absence of water vapour 
makes a difference in the mechanism and kinetics of reaction, 
there will be no difference in the mechanism once the water 
vapour content in the oxidizing medium .exceeds the limit of 
a few parts per million, ITormal variations, therefore, in 
the level of humidity of the reacting gas medium will not be 
a factor in the present experimental work. 
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EPgECT OP SOLIID IMHJBITIES AND CRYSTAL DEFECTS 

Even spectroscopically pure carlon contains impurities 
like boron, titanium, calcium, sodium, silicon, iron, magnesium 
and copper, whicb, can be detected by special means. If these 
are incorporated in the network, the electron densities of 
the neighbouring carbon atoms will be perturbed, so as to 
favour or disfavour attack by oxygen. It, therefore, becomes 
necessary to study the effect of these solids on the mechanism 
and kinetics of oxidation. Traditionally, it has been customary 
to treat the defects in graphite crystals as some solid impuri- 
ties in the structure. This is probably because crystal defects 
and solid impurities have similar effects on oxidation; and 
in some cases it has not been possible to say one from the 
other [23]. Experiments showing a catalytic influence of solid 
impurities do not in general establish how they are bounded to 
the carbon [ 90 ] , However, the fact that gaseous . impurities 
such as oxygen, hydrogen and nitrogen a.toms can be directly 
linked to the carbon atoms whereas solid impurities are indirectly 
linked, probably through acidic groups like -C-OH, is well 
established. It is also generally agreed that most of the 
common solid impurities catalyse oxidation, and substances 
like phosjjhates, borates, pyrolytic carbon, phenolic resins 
and metal carbides inhibit the reaction. An exaustive. list of 
references on catalytic activities is given by Walker [2] 
and Ubbelohde [ 23 ]. Most of these investigators have taken a 
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a coiiiiaon approach of correlating catalytic activity with 
"nohility” or ionisation potential of the elements, obviously 
vvithout much success. Synergistic catalytic effects are invar' 
ably observed v/hen combinations of a number of elements is 
tried. 

The main concern of this review is for the uniformity 
of samples with respect to composition, electronic state of 
impurities, distribution of crystal defects and crystallite 
orientation, both perpendicular and parallel, Fortunately, in 
macro -oxidation, the number of nonuniformities is so larga and 
randomly distributed that on an average the reacting surface 

particles are, two para- 
meters which are badly 
mixed up with the effect 
of crystal defects, impuri 
ties and the irregular 
variation in the effective 
ness factors, Armington 
[80] selected carbon 
blacks of a 7;ide range of 
particle size and re- 
presented the rate data 
as shown in figure 2-6, 

He suggests that crystallite edges serve as zones of high 
resistance to electron flow. The number of edge carbon atoms 


retains its uniformity. 



Seactivity & Specific 
Surface Area 

FIGURE 2-6 
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per ijinit volume, and hence the resistance to the spread of 
electronic effect of an impurity, is more in the case of a 
smaller crystallite size than in a larger one. 

Recently H, Harker and associates [Slj 92] have used 
E.S.R, techniques to study the effect of impurities that 
are added to the graphite net work hefore and after irradiation, 
Hennig’s etchpit technique [60a] has also been used to isolate 
the effect of a defect from that of an impurity. As has already 
been mentioned [11 ], microcinemato graphic techniques have been 
utilized for similar purposes. With the advent of more power- 
ful microscopes, it is hoped that the effect of impurities v;ill 
he better understood, thus paving a way for a clearer under- 
standing of the mechanism, 

ROLB OR MAES TEARSPORT IN &EAPHI1E OXEMTIOR 

The steps involved in graphite oxidation are essentially 
the same as the ones mentioned earlier for general hetero- 
geneous reactions, and. thus the reaction can take place in 
various control regimes. The apparent values of the order of 
reaction, activation energy and the effectiveness factor for 
different zones of temperature are given in Table 2-2, 

Peterson [46], for the first time, using an approach 
similar to the one of Thiele for catalytic reactions, has 
mathematically, treated a single cylindrical pore and a porous 
solid initially containing uniform cylindrical pores with 
random intersections. In essentially a similar manner. 
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TABLE 2-2 


Order of Effective- R-elation 

Reaction ness factor between Remarks 

E^ and E, 

U, Ti 


Zone I 

(lo¥/ temp- ZERO ONE 

eratnre) 


Zone II' 

, (Intermediate >-HALE >HALE 
temperature) 


Zone III - ■ 

(Higli temp- ORE AIMOST 
erature ) ZERO 


The diffusion res is- - 
E = 1. tance and hence acti- 
vation energy for 
diffusion is ZERO 


E^ = E./2 Relation -has short- 
comings hecause of 
missing kinetic & 
diffusion parameters 


E^ ^ E, Both diffusion and 

^ reaction have low 

activation energies 

and hence ma,gnitude 

of E„ is small 
a 


E,B. Eg^ and E^ are apparent and true activation energies 

Walker et.al, [2] gave a formal definition of the "Thiele 
utilization factor" and have derived an equation to predict 
the criteria for the occurrence of a gas-carbon reaction in 
Zone II. Equations for a plane sphero and a cylinder are 
developed on the lines given by Aris'[93] . The treatment 
assumes a first order reaction on the surface of a cylinder, of 
infinite length. A suasisteady state is assumed-. The vari- 
ation of concentration of gas is considered within the solid 
and not in, the gas phase, Yolumetric specific rate constant 
Ky and the effective, diffusivity are assumed constant 

throughout,' the rod. The equation for simultaneous diffusion 
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and chemical reaction in cylindrical co-ordinate system is 
solved for the rate of loss of vjeight of the cylinder. The 
’’Thiele "atilization factor” /rij defined as the ratio of the 
actual rate oi reaction’ to that which v/ould occur if the 
reacting gas concentration were ■uniform throughout the material, 
is expressed in terms of a dimensionless ratio 9 = 
which is an empirical index assuming different values in 
different zones. Results are summarized in Table 2-3, 

TJIBLE 2-3 


Criteria for the prediction of gas carbon reactions entering 
Zone II for samples of various geometry 






Rate of reaction 

: Rate of reac- 



T] xor 

per unit area of; 

tion per unit 

Geometric 


9 > 

A<PII 

exterior surface 

area of exter- 

shape 

9 j 


•for uniform gas 

lor surface 



[ 

concentration 

when 




9ll 

[throughout the 
j sample. 

,9>9a^ 

T 

Plane : 

/■ 

: : 



- - 

thickness Rj 

Cylinders 
radius R. 

Sphere s 
radius R 



2 

w 

|~=Cj^ v^/^ef f 


2/9 

4 



R V^-y/^ef f 

5/ ^ 

6 

^ - — — - — 



Rote that the value of for start of zone II is quite 
arbitrary. In fact, the changeover from zone I to zone II 
can cover a considerable range of temperature. This difficulty 
is not experienced when the reaction passes from zone II to 
zone III, because the overall rate does not vary appreciably due 
to small changes in diffusivitiesi 
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altering this basic background, a large number 
of inTestigators have studied graphite oxidation under varying 
velocities, gas temperature, osygen concentration, etc. 

Many a time air at varying velocities, is admitted over the 
sample and the rate of weightloss is determined. Extrapolation 
in the region of high velocity gives the rate data which is 
free from diffusion effects. After considering the heat transfer 
analogy, and after eliminating the effects of variations in 
viscosity, density, diffusivity and reactivity, many investi- 
gators have attempted to correlate the variation of reaction 
rate with linear gas flow rate as given by the equation 

= a where a is a constant, ■ V is the velocity of gas, 

R is the characteristic dimension of the sample, and n is the 
exponent which can take values from 0.20 to 0.50, Parker and 
Hottel [94], reating birush carbon with air at 1227^0 have 'found 
that the rate varies with 0,37 power of velocity, Mayers [951 f 
using 40- by 60- mesh coke in 1" high beds, obtains a value of 
0,5 for the exponent | Chukhanov and Karsavina [9^ in their high 
velocity experiments using beds of particles 3 by 5.5 mm in 
diameter foimd a value of 0,4. Kuchta and Co-workers [ 97] using 
carbonrods report an exponent of 0,45, Pay [78] using carbon 
and graphite rods reports a value of 0.5, and lu et.al, [98] 
a value of 0,49* More and more experiments were done because of 
the differences in the value of exponent n. However, it was 
Graham [64] v/ho studied the steam-carbon reaction under high 
velocity conditions and clearly pointed out that an exponent 
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less tlian 0.5 tueans that the reaction is not in zone III and 
an exponent of 0,5 indicates that the reaction is well in the 
high temperature zone. It also goes without saying that Under 
identical reaction conditions (i.e, for a fixed composition 
and velocity of the gas stream at a fixed temperature) the 
rate of reaction in zone III is independent of the type of 
carbon reacted. Some reservation can be expressed regarding 
the nature of carbon effecting the carbon-^honoxide-carborf -dioxide 
ratio lea.ving the surface and hence the reaction rate per unit 
of oxygen, diffusing to the surface. Even though sufficient 
data is not available, Day [78] reports that the reaction rates 
of petroleum coke, graphitized lampblack, and graphitized 
anthracite rods agree within 12^ at a temperature of 1827°C 
and at a constant gas velocity, 

Por reaction at the same temperature, P, I, Walker [20 has 
theoretically predicted the relative rates of the different 
gas-carbon reactions, in the high temperature zone, when 
using a sample of fixed dimensions, a constant linear gas 
velocity and a fixed concentration of reacting gas in the main 
stream, (Refer Table 2-5), logically the relative rates 
should depend on 

(i) viscosity and density of gas present in the stagnant 
film, i,G, the parameters occurring in the expression for film 
thickness, 

(ii) free diffusivity through the stagnant film. 
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(iii) num’ber of moles of carbon reacting witli each 
Qole of the reacting species, . 

■The reaction rate constant does not come into the picture 
because oi tne diffusion controlled regime. Surprisingly, 
very few investigators [ 52 ] have tried to exploit this idea 
to determine one of the above four parameters when the 
remaining throe are 

Occassionally, the change in physical structure of 
graphite, in terms of densr'-r and area profiles, after degrees 
of burn off a,t different temperatures can aid in the unde 2 >- 
stonding of gas-carbon reactions, especially so in deciding, 
the zone of the reaction. This is well illustrated by ¥alker[2] , 
Sorae investigators [60b, 29 J have tried to correlate the poro- 
sity and the poresize distribution with the percentage degree 
of burnoff of the sample, Figures 2-7 and 2-8 indicate interest- 
ing results of Walker for burning of hollow cylinders of graph! t( 
Hov/over* a satisfactory explanation of the extremum points in 
the curves is not available. Similar multiple extremes seen 
in the investigations of Griffiths [21] and Thomas [29] are 
also not well explained. 

Most of the recent investigations in the field of mass 
transfer with chemical reaction of graphite, have been under- 
taken with one of the following obijectives in mind; 

(a) To evaluate the properties of graphite as a 
material of construction for rocket nozzles. 
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(b) To improve the performance of graphite as a 
moderator in nuclear reactors, 

(c) To estimate the time dependent effectiveness factor 
for none ata.ly tic reactions, 

P, Havjtin and R,A, Huber [ 52] have studied experimental oxi- 
dation of a tubular specimen of nuclear grade graphite in 
v^hich oxygen access is confined to the boro of the tube and 
the mathematical model developed has been used to estimate 
chemical reactivity, utilization efficiency and the effective 
diffusion coefficient through the porous ■walls. In fact, this 
work- looks like. the starting point for the present work, II,G, 
Schweitzer, R,M. Singer and Co-authors [ 51, 99] have studied 

r 

the siQultane.:ais oxidation and heat transfer in graphite 
channels with an obvious objective of predicting the stable 
performance of air cooled BHI graphites in nuclear reactors. 
They have covered the 'experimental and theoretical aspects of 
parameters like temperature, flow rate, channel diameter and 
chemical roactivity. Two important conclusions arrived at 
are relevant to the present work, firstly, the most rapid 
temperature Instability is concluded ‘to be due to the secondary 
gas phase . reaction of carbonmonoxide with oxygen. Secondly, 
no diameter effect other than that accounted for by the heat 
transfer coefficient was observed implying thereby the absence 
of radia,! concentration gradients in the system, Maberry 
et.al , [28] have tried to evaluate the performance of a graphite 
nozzle subjected to internal errosion and corrosion. Straight- 
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forward dimensional analysis is carried oni: "by nsing several 
parameters like nozzle diameter, length, gasflow rates, fuel 
to oxygen ratio, mach number, temperature, tliroat to mouth 
pressure ratio, and the diameter ratios. The exponents for 
various dimensionless groups are estimated. As late as in 
1970 [1003, aolation performance of various types of graphites 
have been tested with a view to assess the potentialities 
of these materials for the use in dynamic aerospace environ~ 
ments. Hemispherical samples during the course of diffusion 
coiitrolled reaction, turn out into conical shapes, but unfor- 
tunately no mathematical model is proposed. The models given 
by Sc ala [1OI, 102] do not agree with the experimental results. 

Blaahs [ IOO] attributes this to a pxobeble reduction in 

reaction rate at high temperatures than the ones predicted 

by the extrcipolatod Arrhenius plot, G-.L, Tingey et,al, [ 30] 

have studied the thermal and radio lytic reactions of coolant 

impurities with graphite in a high temperature gas cooled 

reactor. The Bessel function solution of an equation for 

diffusion-limited reaction is numerically solved on the 

computer. Even though the temperature variations along the 

length and radius are considered, variations in v/all profiles 

have not been determined owing to an extremely slow rate of , 

gasification of the moderator, Ralph ITorman [103] has treated 

the problem of mass loss as a function of time and position, 

but the sanple geometry is a flat plate and that too with an aasumptioit; 

of one dimensional slug flow and a quasisteady state condition. 
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3. SIATERIADS MD METHODS 


Before going info the details of the experimental 
proced'are, it may be worth mentioning the specifications 
01 the materials vsed, the list of which is given below: 

(1) GEAPHITB: Metallurgical grade graphite which had been 
used as a melting electrode in a carbon-steel manufacturing 
furnace for over, an year, was used as the starting material 
for preparing the samples. The apparent density measured 
by weighing an accurately machined cylinder was found to 
vary from 1,600 grams/cc to 1,605 grams/cc. Depending on 
the temperature of graphatization, the theoretical density 
of nonporous graphite is known to vary between 2,26 gram/cc 
and 2.10 grams/cc. Based on this, the variation in the 
porosity of the solid was estimated to be in the range of. 
0,29 to 0 , 23 , The ash content determined by the standard 
A,S,T,M, method [104] wr.as found to be less than five hundred 
parts per million. 

(2) IflETHAHOL: EDH analar grade solvent ^ 

( 3 ) PLASTER OP PARIS: No particular specifications, 

( 4 ) SILVER GOOP: A patented paste of a silver compound, 

(5) NICKEL POILSs Obtained from 99.99°/ pure nickel sheets,. 

(6) SILVER WIRE; Made of 99.99/ pure aetalic silver, 

( 7 ) SILICA: Pine ly powdered form with no particular speci- 
fications. 
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(8) !FELSPARs Finely powdered form with no particular speci- 
fications. 

(9) GHINA CLAYo Finely powdered form with no particular 
specifications, 

(10) CALCITO OXIDE; Finely po¥/dered form with no particular 
specifications., 

(11) SOROSIL CSJIMT; A patented paste of a special quality 
cement dissolved prohably in a mineral oil, 

(12) SEDZOIC ACID; Fine flakes of B.F, grade. 

(13) SODILRI HIDED FIDE; Solutions of various strength, 

(14) REFRACTORY BRICKS; Fire proof alumina bricks capable 
of withstanding temperatures upto 1800°C, 

(15) GOliPRSSSED iilR; As available from a compressor or a 
cylinder, 

(16) ’OrYGEN; 959^ purity gas. filled in a cylinder, 

(17) KITROGEK ; 98^ purity gas stored in a cylinder, 

(18) DIQUID KITROGEU; As available from the liquid air plan. t. 

The entire experimental work can be divided into 
three parts; 

(a) Oxidation of graphite in a single pore, 

(b) Reaction of benzoic acid and caustic soda in a 
pore, 

(c) Determination of the first order reaction rate 
constant. 
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ErPERlMTAL PROCEDURE FOR (fi) 


^ Qj — §.-iyiLSi£!^ • G-raphite blocks available in irregular 
sizes and siiapes were macliined on a precision latbe to obtain 
o, number of be llow cylinders tbe dimensions of wbicb are 
indicated ir ■"'ible 3-1. 


TABLE 3-1 


Dimension in Cms, 

Set 1 

Set 2 

Set 3 

Length 

o 

t 9 9 

OJ 

4.0 

6.0 

Inteina.1 diameter 

0.1 

0.2 

0.4 

External diameter 

0.5 

1.0 

2,0 


Some of tbe guiding principles and constraints wbicb 
restricted tbe maximum and minimum dimensions of length, 
inteinal diameter and external diameter, are enumerated 
below, Tbe limiting factor for the maximum length was tbe 
operating length of the available precision drills. Holes 
were drilled from both tbe ends, and the minor eccentricity 
resulting therefrom was corrected by sliding and rotating 

It 

the cylinder over a brass rod of accurate dimensions. In 
fact, drilling was started with a bore- of a smaller size 
and was gradually increased to a diameter slightly less 
than the desired value. The minimum length was restricted 
by the desirability of having an appreciable gas phase 
concentration gradient in the axial direction, ibr the same 
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reasons, tlae aaximuia internal diameter had to be restricted 
to 0,4 cm. As regards the minimum internal diameter, the 
limiting factors were an ‘unduly long reaction time, a long 
machine time for computation, and the Tibrational stresses 
causing the "breakage of small diameter drills, There were 
no serious limitations experienced in the choice of the 
external diameter of the cylinder. While a larger external 
diameter was suitable for the purpose of coating and material 
handling, a smaller one was preferable for a short time of 
reaction and computation, 

PBETARATIOH OF BEACTII-IG SURFACE ; To obtain uniform surface 
properties of graphite, samples were v/ashed with methyl 
alcohol [80] « Machined graphite samples were placed in a 
conical flask and acid washed methyl alcohol was poured 
into it. The flask was corked and placed in a constant 
temperature bath at 40°C for about ten minutes. The samples 
were removed and dried in an" oven at 100°C for two hours; 
they v/ere now ready for a protective coating. At this stage, 
the weight of all the samples was recorded. The time lag 
between the surface preparation amd the mounting of the 
samples in the furnace was not an important fa,ctor affecting 
the initial rate of reaction, 

COATIMGr OF EXTEBIAL SUBPAGE ; The technical problems involved 
in the insulation of the external surface against chemical 
reaction ■were many and varied. The various methods tried for 
obtaining a stable protective coating are listed below; 
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(i) A pastie of Plaster of Piaris was coated on tiie 
excemal surface ard Y;as allowed to "set" in an oven at 
temperatures varying from 80°C to 300°C. The coating thus 
obtained v/as porous and could not act as an airti^t covering. 

As soon as the samples were introduced into the furnace, 
cracj£s in the coating became prominent and the samples 
sta,rted burning from outside, at a rate which was much faster 
than that of the reaction inside the pore. The method, there- 
fore, had to be abandoned. 

(ii) Some silver compounds are known 1^23] to act as 
inhibitors for graphite oxidation reactions, A patented 
material "Silver G-oop" was coated on the external surface 
and the sample was dried in an oven at 200®C for two hours. 

Again, hov/ever, when the sample was introduced into the 
furnace, the coating became unstable and allowed the external 
surface to react; and this method also had to be discarded* 

(iii) Attempts were made to generate a coating of 
silicon carbldo on the external surface so that it remained 
stable and unreacted even 'upto a temperature of i8G0®C, 

Solid cylinders were completely covered by finely powdered 
silica and were introduced into a furnace, Carbide formation 
normally occurs at temperatures above 1800°C but the porous 
structure of graphite began crumbling even at 1650®G; so 

no attempts were made to reach a temperature of 1 800° G for 
obtaining a silicon carbide covering, 

(iv) It was contemplated that a metallic sheathing jj2m could 

^ ^ ! nj yfr 
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act as a protective coating, Tlie graphite cylinder was 
inserted into a stainless steel However j. the co- 

efficient of expansion of metal being higher than that of 
graphite, a tight fitting could not be retained during the 
course of reaction. Moreover,^ softening and distortion were 
added problems. In fact, the softened metal, on entering 
the latix of graphite, catalysed the surface oxidation in- 
stead of inhibiting it. Protection by means of a metallic 
tube hence v/as out of consideration. 

In order to avoid the problems of softening and dis- 
tortion of a thick metal coating, attempts vjeve made to deposit 
very thin layers of nickel and silver vapours on the external 
surface. Graphite cylinders were rigidly mounted on rotatable 
shafts connected by aluminium gearwheels, one of which was 
j^igidly fixed on to the shaft of a micromotor. The system was 
so placed that the tungsten boat of the Yeeco vacuum evaporator 
was centrally located at a point equidistant from the four 
graphite cylinders, which could freely rotate about their 
respective axes. Ihe metal foil was placed on the tungsten 
boat. Mter closing .the belljar of the Yeeco evaporator, the 

two vacuum pumps were started in succession, and v/hen the 

— »6 —8 

pressure was 10 to 10 cm,, of Hg. , an electric . current was 
passed through the metal foil. Then the micromotor was switched 
on so that the graphite cylinders began to rotate, thus 
enabling a 'uniform coating^of the metal vapour on the external 
■surface. Contrary to our expectations, such a layer did not 
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adhere to the surface; this in turn rendered the method 
ineffective, 

(v) hext, ceramic protective layers such as frit and 
glass were tried. Materials were intimately mixed in the 
proportions indicated in Table 3-2. 


TJffiLE 3-2 


Subs tanc e 

Percentage Composition 

V/eight in 

Grams 

SET 1 

SET 2 

SET 1 

SET 2 

Silicon dioxide 

80^0 

80/c 

200,0 

160.0 

felspar 

1 05^0 

1 

25.0 

30.0 

China Clay 


Afo 

22.5 

8.0 

Calcium Oxide 

\io 


2,5 

2.0 

TOTAL 

100/o 

100^0 

250.0 

200,0 


The mixtiire of composition mentioned in set 1 was 
mixed with 300 cc of water and poured into a stone ware ball 
mill which was then closed and placed on rotating bars. 
Grinding wa.s continued at 60 r.p.m, for about twenty four 
hours. The resulting frit was used for coating the external 
surface of the graphite cylinders, which after drying for 
8 hours at room temperature, were introduced into a furnace 
at 600°C, Nitrogen atmosphere was maintained by constant 
flow of nitrogen through the furnace. The temperature was 
allov/ed to rise to 1000° C. in an hour ’s time. However, the 
coating thus obtained had developed cracks in it. The frit 
composition, therefore, was changed to the one mentioned in 
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set 2, but this also did not yield any favourable results, 

(vi) Lastly, a patented paste of "Sorosil Cement" 
was tried. This coating appeared to be stable and was partially 
successful in preventing surface oxidation. When the initial 
adherance was proper, the coating used to serve its purpose 
but on many occasions, due to improper adherance, the coating 
allowed the diffusion- of oxygen into the space between the 
outer surface of the sample and the inner surface of the 
coating. This naturally caused the oxidation on the external 
surface to produce burnt samples of the shape indicated in 
Lig.5-1* In such cases, since the rate of burning on the 
outer surface was much higher than that inside, the results 
for reaction in the pore could not be obtained satisfactorily’’ 
However, when the coating was stable, the sa.mple was 
mounted on to a fire proof brick and placed in a furnace in 
such a way that the axis of the cylinder was perpendicular 
to the direction of flow of the oxidising medium and parallel 
to the horizontal direction. The flow of oxygen and nitrogen 
v/as adjusted in such a way that the carbon dioxide content 
in the exit gases remained less than 2 .$^ The outlet composition 
was checked from time to time by means of a gas chromatograph, 
tAfter a desired interval of time, the samples were removed 
from the fu 2 ?nace and allov;ed to cool. They were then weighed 
with their coatings. By gentle tapping, the protective shells 
were broken., and the. resulting specimens . were cut into a 
number of rings to measure the radii at various lengths and 
time intervals. 
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Only one of the ahove experimental rnns is { 

I 

reported in the chapter on "Calculations and Conclusions" 
for the following reasons? 

(a) Hardly one out of six sorosil cement coatings 
adhered properly and hence the results were neitii'er:' repro- 
duceahhe nor relia.hle, 

(Td) The available furnace - even if operated well 
below its maximum rated temperature ~ could not be used 
continuously for more than 36-48 hours, without running 
the risk of burning the heating coils,. In many cases, meaning 
ful results could be obtained only after a long time of 
reaction. 

(c) While weight loss measurements were somewhat 
dependable, profile measurements were not accijtrate enough 
have • a meaningful comparison with the theoretical values, 

EXEBRIMMIAl PROCEhimE fOR (b) 

In the face of the above difficulties, the possibility 
of some low temperature reactions was explored. With a view 
to obtain a hollow cylinder of an organic material, a brass 
mould indicated in Eigure 3-2 was fabricated, 1 groove was , 
canved oh the base plate (a) to, allow a pushfit of the 
brass pipe (b) $he threaded portion of the central rod (c) 
could be fixed in the central hole of the base plate, A 
layer of silicone oil was applied on the surface of each 
part of the mould. Molten organic substance was poured into 
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tliG annular space of the mould assembly, which was allov>?ed 
to cool and solidify at room temperature. The cylinder 
thus prepared v/as then removed by slow rotary motion of the 
brass tube and rod. 

The above brass mould was used to prepare hollow 
cylinders of benzoic acid for on expex’imenta,! illustration 
of low temperature solid-liquid reaction in a pore, How’ever, 
the homogeneity of the sample could not be ensured because 
of nbniHiiform changes in the chemical structure of benzoic 
acid during the process of moulding, and because of an uneven 
distribution of air gaps in the solid phase of benzoic acid. 

The reaction of benzoic acid with caustic soda was not consi- 
dered for the main stream of experimentation because it was 
felt that this combination of low diffxisivity (10*"^ to lO^^cth/ 
sec) and high rate constant, may not give a good spread of 
results, Iforeover this diffusivity brings in very low value 
of the coofficient in the quasilinear parabolic differential 
equation, introducing serious instability in numerical solution. 
However, 

(a) selection of a suitable low temperature system, 

(b) investigating suitable ranges of diffusivity and rate 
constant, 

(c) exploring Dankweirt’s criteria for instantaneous reactions, 

(d) perfecting a method of moulding solid cylinders in inert 

atmospheres, ; 

(e) trying pressure pelletization techni< 3 .tie for solid moulding. 
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(f) choosing a proper range of temperature to realize "both 
diffusion and reaction rate controlled regimes, 

(g) selection of a proper coating for the external surface, 
and 

(h) exploring the possibility of a homogeneous catalytic 
reaction, can form a potential field for future research 
investigations in the area of noncatalytic reactions in 
porous solids in general a single pore in particular. 

EXPEBIllBITM PRQGEDITBE FOR (c) 

In spite of the difficulties experienced in a high 
temperature operation, and may be because of the limitations 
of a liquid-solid system which are already indicated in the 
introductory note, persuit for standardization of a suitable 
method for the determination of the reaction rate constant 
was continued. Traditionally a flow system is preferred for 
the purpose, but in the present situation a static system 
has been selected with a viev; to simulate the conditions 
prevailing in a single pore. 

Machined solid cylinders of graphite were cut to 
obtain a. number of cylindrical tablets. Each tablet was 
rubbed on the surface of a f'irte sand paper to reduce its 
thickness to less than 2 millimeters so that the reactions 
on the curved surface could be neglected. The reacting 
surfaces of all these samples v;ere again treated with acid 
washed methyl alcohol as described before. The dried samples 
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v/ere weighed the cross-sectional area of each tablet was 
measured hy means of a calliper. The samples thus prepared 
were now ready for reaction. 

The furnace associated with the T.G^A. (Thermo Gravi- 
metric Analyser) was switched on for obtaining the desired 
temperature. The weight loss axis was adjusted to cover a 
span of 500 mg. Depending on the temperature, the time axis 
was adjusted to cover a range of half hour to three hours. 

Cooling water was circulated, Finely powdered silica was 
filled in a porcelain crucible; the sample Y/as placed on the 
pov/der and gently pushed inside so that only one circular 
face was exposed to the atmosphere.* The crucible was then 
placed in a vertical glass casing which formed an integral 
part of the T.G.A, Counter weights on the (Quartz spring V7ere 
adjusted to obtain the desired position of the pointer on the 
graph to indicate the initial weight of the sample. As soon 
as the furnace attained the desired temperature,, the gl^s 
casing was lowered into the furnace and the temperature over 
the sample was noted at half minute intervals. Vftien a temperature 
of 600®C V7as reached, time and weight recorders were switched 
on and temperature was noted at regular intervals of time 
until the sample attained the temperature of the furnace^ 
which was set at a particular value to be taken as the reaction 
temperature. After permissible levels for weight loss and 
time were reached, the sample was removed and checked for 
possible reaction on its other face which was totally covered 
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"by silicon dioxide* llie occnrance of reaction on the covered 
surface could he judged hy visual observations, and -when the 
sample was found to have been excessively attacked, that 
particular experiment was repeated. The weight loss curves 
for temperatures between 600°C. and 850°C are shown in figures 8—3? 3**^ 
^At temperatures above 850°C, weight loss recordings 
on T,G,A, could not be done because of the prolonged severe 
noniso thermal conditions that prevailed in the initial stages 
of reaction. On the other hand, if the samples were intro- 
duced directly into the furnace without weight loss recording 
devices, the time for which the nonisothermal conditions existed, 
i,e, the unsteady state time, was much less, for a hi^er 
temperature zone, therefore, this kind of furnace was used, 
and weightless measurements were mannually done at discrete 
intervals of time. 

Blind holes, each of one and a half inches depth, v/ere 
drilled into pieces of refractory bricks;, and silica powder 
was filled into them so as to form a layer of one half inch 
thickness on hie bottom surface. An accurately weighed sample 
was placed on the powder and was gently pushed inside-, so that 
only one face was exposed to the atmosphere. Eight to twelve 
identical samples of this kind were simultaneously introduced 
Into a furnace, the temperature of which was set at the desired 
value. As before, compressed air was passed continuously 
through the furnace so that the carbon monoxide and carbon 
/dioxide content in the exit gases did not exceed 2?^,- Exit gas 
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composition was claecked by means of a cliromatograpli from 
time to time* After every 15-20 minates of reaction, one 
sample was drawn out of the furnace; wbb cooled and then 
v/eighed^ l/'eight loss curves drawn from this data are shown ' 
in figures 3-.5,3-.6(p« .82 & 81). More reliable data regarding' 
numerical values of rate constants and limiting slopes was 
obtained at low temperatures than at high temperatures because 
tlie parameters temperature and weight were continuously 
recorded in the former case. 





■4-. MATimmiCKL MOIDELLIIG 


There are two aspects of the mathematical formulation 
of the problem under consideration, the first One dealing 
with the determination of the first order rate constant for 
a surface reaction and the second one with the main problem 
of a heterogeneous fluid“Solid resicticn in a single pore. 

Before going into the details of those aspects, it may be 
worthwhile enumerating the assumptions involved in the 
mathematical analysis of the system, 

(l) ISOTHEEM-AIi COUBITIOh i There are four possible situations 
wherein nonisothermal behaviour is likely to occurt : 

a) At the beginning of the reaction, the temperature ■ I 
of the, sample increases from room temperature to 

the temperature of the furnace, 

b) Similarly at the end of the reaction the sample 
cools, from furnace temperature to room temperature, 

c) Dissipation of the heat of reaction nay result in 
radial temperature gradients within the solid. 

d) Oxygen concentration gradients may cause temperature 
gradients in the gas phase. 

A q.uantitative or a qualitative analysis of each of these 
nonisothermal situations will now be made to conclude their 
insignificance in the system. 

To begin with, the time required to raise the temperature 
of the graphite cylinder from room temperature to furnace 
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temperature is estimated* The most severe conditiorfe that 
are likely to occur are assumed. Even though the effective 
thermal conductivity of -graphite increases with temperature, 
its minimum value is taken for calculations. Also, the 
maximum numerical values for reaction temperature, specific 
heat and cylinder radius are assumed.. The heating time 
can he estimated as followss 

Effective ' thermal conductivity for graphite K = 0,5 
cal./cm.sec,°C [105] ‘ ® 


Radius of solid graphite cylinder R = 2,0 cm. 

Room temperature = 20 °C 

Eurna.ce temperature = 1550®C 

Length of the solid graphite cylinder 2H = 6 cm. 

Density of the graphite used = 1,6 gm/cc 

Specific heat of graphite = 2,5/l2 cal/gm,°C 

ps 

Thermal diffusivity *^ps^ 



2 

= 0,9 cm / sec. 


With the help of a standard Gurney Lurie chart for unsteady 
state heat conduction, the time required to heat the sample 
from room, temperature to 9 9. -95^ of the furnace temperature is 
estimated to he fonir seconds . This being true of a solid 
cylinder, the time required for a hollow cylinder will naturally 
he still less, Ro separate calculations need he carried out 
for the cooling period because it is not likely to exceed the 
heating period, , ^ ^ 
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In the actual experiment,, however, the heating and 
cooling periods varied from 80 to 1-35 seconds. This difference 
between theory and practice might have been due to the extra 
resistance to heat flow offered by the protective surface 
coating over the gra,phite cylinder, in important fact that 
is revealed by the above estimates is that the heating or 
cooling period of a few minutes is much less than the reaction 
time of a fe?; hours. Besides, the temperature does not remain 
at the maximum level throughout the unsteady state, and hence 
the error introduced by neglecting the unsteady state time is 
not as severe as it would have been if the temperature had 
remained at its maximum level throughout. The only drawback ' 
of this assumption is that? the severity of error increases with 
temperature-. At high temperatures, the total time of reaction 
?/ill be less and the period of cinsteady state for initial 
heating and final cooling will be more so that the percentage 
errop will be higher. Taking into consideration all these 
aspects of unsteady state, the initial heating period is 
included in the reaction time while the final cooling period 
is excluded from the same, . ■ 

Coming to noniso thermal behaviour due to the heat of 
combustion of graphite., the maximum possible rise in temperature 
of solid can be calculated as follov^ss Again, sissuming. 

Equivalent thermal conductivity 'of graphite,!; = 0,3 cal/ 

sec, cm,°C, ® 

Maximum heat of combustion, 0 H = 7831’ cal/gm of graphite. 



. iPirst order rate constant, Ic = 10,0 cci/seG,> 

Mariaum concentration of oxygen in air =.21^, 
density of air = 0,20 x 4-#464 x'^IO"^ gm- moles/cc, 
Externa,! radius of the cylinder Ee = 1,0 cm., 

'Internal radius of the cylinder Rj^-= 0,2 cm., 
and denoting the external and internal surface area of the 
cylinder by k and A. respectively, the standard steady state 
heat transfer equation for a cylinder is given .by 


' A ^ -A. 

(12)(2 7cR, )( OH h G. ) ^ = K ■ 

1 1 A max e (A 

T 1,^ / 


X 


4 I 

. 


, , , 4""1 


ln( 


Therefore, the maximum possible temperature rise = 1.95°C, 
Including the resistance offered by the protective coating, O ® 
is not likely to exceed the value of 4°-5°C, 

In these calculations, heat transfer in the axial 
direction is not taken into consideration. This error, however, 
is on the safer side because if axial transfer is considered, 
the rate of heat dissipation will be higher. Moreover, the 
above calculations being- valid for most severe conditions, the 
error in an actual experiment is likely to be less severe than 
that estima.ted by equation 4-1, The internal surface area 
increases during the course of reaction,, but at the Same time, 
the thickness of the cylinder and oxygen concentration decrease, 
so that any adverse effect of the increased reacting surface 
area is nullified. Thus, from theoretical considerations, 
it is seen that the noniso thermal condition created by the rise 
in temperature is insignificant. 
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In tiiie actual experiment., lioweTcr, the thermo couple 
placed in the inner cavity did not register any a.ppreciahle 
change in temperature during the course of reaction, This 
was because the effect of sluggish heat conduction through 
the graphite wall has a tendency to counter' the nonisothermal 
effect of heat of comhiistion. Theoretically, it takes infinite 
time to raise the temperature of the sample from 99.9?^ io 100?S 
of the furnace value. 

Lastly, the possibility of the oxistance of an axial 
temperature gradient resulting from the corresponding axial 
concentration gradient is considered. However, in the present 
sitiiation, the teraperatu.re gradient is unlikely to occur 

bejCause the area for heat transfer is much more than that for 

> • _ ■ • 

I 

mass transfer, as indicated in figure 4-1, In other vrords, 
the rate of heat dissipation is much higher than that of m^s 
transfer, so that the temperature gradient, if established,. 

Vi/ill tend to vanish much faster. Even in the actual experimentj 
no change in temperature was observed when the thermocouple 
junction was shifted in the axial direction. Therefore, - the 
assumption of isothermal condition in genera.1 is not likely 
to introduce any significant errors in the mathematical analysis 
of non-catalytic reaction in a pore, 

(2) GOKDITIOH OE COHSTMT DIEECTSIVITY s The Isothermal condition 
having been justified, the variation of L need be considered 
only with respect to the concentration of gas. At room 
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p 

■bemperature , the diffusivity of ozygen in air is 0.2 cm /seo 

p 

and tii 8 .t in pure c art on dioxide is 0.18 om /sec [l!i06] . While 
these iiumericod. vahues may partially prove the validity of I 

the assuraptioii of constant diffusivity, it is extremely diffi- 
cult to provide an eJialytical .justification, I'liia is because i 
ilickss la?/ of diffusion in itself is an over simplified 
ma,thematical description of the actual process of molecular 
diffusion, A constant numerical value for the diffusion co- 
efficient is more often a mathematical necessity than a physical i 
reality. In such situations the tradition established by 
various investigators becomes a good justification for the ; 

assumption, 

( 3 ) ABSBilOS Qg QOhYEGIIOI'T ; Wlien the axis of the cylinder 
lies in the horizoP-tal position, there are two ways in which i 

convection currents could be set up: | 

(s-) Gravity ; Usuallrf this effect is negligible over 
small height se and since there is a radial symmetry of i 

temperature and ooncentration, natural convection is non..^existent 
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possibility of convect.iori in tihe axial direction does 
not exist becanse an axial concentration gradient can cause 
only molecular diffusion., 

(b) glow Past Oylindcr s This can set up convective 
currents near the mouth of the pore. Employing low air 
velocities miiiiniises convection. Moreover, ceramic protective 
coatings on the two sides of the cylinder act as dampers for 
air currents. Computation of a momentum entry region for 
flow past a hollow cylinder is rather an involved procedure, 
and is not attempted here.. It may, however, be noted that as 
the cylinder shrinks in its length, the effects of convection 
will be reduced considerably because of the increased distance 
between the pore mouth and the source of disturbance, 

(4) PRODUCT OE OXIDATIOII IS GOp S Sufficient evidence 
for this has a,Iready been given in the revie?-/ of literature. 

Even though the assumption looks simple, it considerably 
reduces mathematical complexity, 

(5) FIRST ORDER REACTION ; That the surface reaction 
between graphite and oxygen is of first order with respect to 
oxygen pressure, has already been concluded in the review 
discussion, ilny deviation from this,, if significant, will be 
discussed in the concluding chapter, 

(6) SURFACE IREA FOR PwEAGTIOD ; During the course of 
reaction, a large number of closed pores open up and there is 

a considerable change in. the' surface area available for reaction. 
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The gectnetric surface ^ea, in fact, v/ill le a negligible 
part of the total reacting surface area. Pirst order rate 
constants are, howeTer, defined on the basis of geometrical 
surface area throughotit the reaction. Various inyestigators[ 2] 
have found this approach satisfactory, and in the present 
analysis also, a similar procedure is followed. Such an 
assumption is well justified by the fact that the limiting 
weight loss curve for each temperature is more or less a 
straight line. Invariably, the period for which the weight ■ • 
loss curves are nonlinear is much less than the total time of 
reaction. I’herefore, it becomes convenient to define rate 
constants on the basis of the linear portion of the weight -loss 
curves and the geometric surface area of the sampled 

(7) RihlhL HEfPgSIOh ; Almost all the authors have 
ignored radial diffusion in the analysis -of a chemical reaction 
in a pore, either in soi implicit manner or by mahing a quali- 
tative reference to the negligible transverse gradients in 
a pore of high length to diameter ratio. However, .no qhanti— . 
tative analysis is reported regarding transverse diffusion ■ flux 
in a reacting pore, A modest beginning in this .direction can. ■ 
be made by starting v/ith a two dimensional unsteady state 
diffusion equation in the cylindrical coordinate system. Let 
the axial distance be denoted by 8, internal radius by ' 

the diffusivity of oxygen in air by D., the time by t, the • 
length of the hollow cylinder by 2H, oxygen concentration by 
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the initial oxygen concentration hy 
Defining a dimensionless axial distance* 
a dimensionless concentration 


X = 

® = '’a/® AO, 


a dimensionless radial distatice 
the two dimensional irasteady state diffusion 


r = R/R^, 
equation may be 


v-/ritten as 


0^1 h , 1 rC 

_2 r ^ 




■t:? r 


D at 




For a quasi steady/ state the equation 4-2 hecoines 




dx' 


-j. 

2 r 


1 .3 , a c 


=0 


« . « 4*^5 


If this pseLido steady state is considred. after a time interval 
the internal wall radius R will, in general, be a function 
of both X and t^ , say R = g(x, t^), The boundary conditions 


\vill then be 

C(o , r) = 1 ’ . . .4-4 

-D (D,r) = 0 (where I = H/R. ) ..«4-5 

■s 2; _ X 

““ (x,0) = 0 ' ...4-6 

“D -~4x, g(x,t,,)J =_k R^C-f x,g(x,t^)j-_ ...4-7 


Ihe functional relation g(x,t^) being unlcnoTO, the boundary 
condition 4-7 cannot be used to integrate the equation 4-3. 
As a first approximation, if a linear wall profile i, e. a 


*In the subsequent- ana lysis, dimensionless axial distance 
is def ined as x=.2/H. This change in notation will be 
referred to as and when it appears. 
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conical pore surface is assumed tlie boundary condition 4~7 
becomes 

“D (x, ax + b) = k R.C(x,ax + b) ^,■.4-8 

wliere a and b are constants for tlie sla.nt beiglit of tke cone. 
Even with the simplified boundary condition 4-8, an analytical 
solution of equation 4-3 does not seem possible,- However, 
if the interna,! ?;all radiiis is independent of axial distance, 
(which is tme in the early stages of reaction), the relation 
4-8 becomes 


Ir ^ ^ Rj_C(x,l) 


...4-9 


By using finite Hankel transform H-echnique, the solution 

in terms of standard Bessel functions, to equation 4-3 with 

boundary conditions 4-4, 4-5, 4-6 and 4-9 may bo obtained as 
00 


2 


a,; 


jpai) J„(o,_r) 


cosh a^(L-x) 


r r N , t 2/ cosh a.L 

[P (a.)+J,,(a.)] 1 


where a^J,j(a,.) = 




JoCaj.) 




and 


I 2kIL 
m = h'xj ■ 


...4-10 

,,.4-11. 

...4-12 


ITote that the dimensionless parameter "m*' represents the 
standard Thiele modulus for noncatalytic reactions. 

To determine the tvro dimensional effectiveness factor, 
it becomes necessary to integrate the gasphase concentration 
both in radial and a,xial directions. Since the pore mouth 
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flux gives the actiis.1 rate of react ion over half the length 
of the pore, integration need he done only in the radial 
direction. Thus average concentration is given hy 


/2 


t:rdrC ,1 

^ 2 I Crdr = 4-) ta, [1+ ■ g 1 ) 

J27Trdr 't a •• (^i) vf 


-1 


cosh a^(Ii~x) 
cosh 1 


...4-13 


Therefore the two dimensional effectiveness factor hocomes 

-1 


-D 


^2 = 


..0 


2^ R.Hk C . 

1 Ao 


-p 

21 


SL.'\“ 

-r3/-- 
h a 




tanh(a^l) 

. . .4-14 


hom^ the error introduced hy neglecting the radial diffusion 
may he estimated hy comparing this £2 with , the one 
dimensional effectiveness factor, i'or this, % useful approxi- 
mation to ^2 values of m/l can he obtained hy 

using the relations 

a^ ^ Cv»-* ^ ^ C ...4— lp 

and a 2 — ^3.8 as g yO ,,,4-1 d 


The equation 4-14 thus redoes to 


.,v tanh ni 
m 


1 


siranf)" + 0 (E)^ 






(1 + 


2^ 1(21)2 Q^ni)3j 

sinh 2w ^1'^ ^1^ 


...4-17 

...4-18 


As expected, Sg become the same as m/l — ^0, and 

from the error term in equation 4-18 it is obvious that 
^2 is always less than This is also in accordance 

with the physical reality that two dimensional diffusion 
reduces the effectiveness of resetnA-n' 
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Is lias Tieen aptly pointed out by Aris [107], simpli- 
fications of the type 4-17 are not valid, "because in the 
present context, typical numerical values of m/l are of the 
order of ten. At least in principle, radial diffusion cannot 
he neglected even for a simple cylindrical geometry. Now, 
what happens when the pore wall profile is not cylindrical is 
eiiy'body's guess, A cautious approach is necessary for any 
intuitive conjecture. 

Inspite of these clestr indications of importance of 
radial diffusion, a large number of investigators have made 
use of the one dimensional effectiveness factor v^ith surprising 
accuracy, lor mathematical simplicity, the present investi- 
gation also assumes no radial diffusion, The possible advantages 
and errors arising out of this assumption in the analysis will 
be discussed in the concluding chapter, 

(8) S.AI'iPhh HETEROGENEITY g As has been pointed out in 
the revieviT, a number of parameters like porosity, po resize 
distribution, nature of bonding between impurity and carbon 
atoms, degree of graphatization, distribution of active 
and defects can and do cause heterogeneity in the sample. 

However, the number of parameters and the number of sites for 
each parameter are so large that the problem of heterogeneity ■ 
acquires a statistical nature and the effect evens out on 
macro scale. Instances to this effect have already been cited 
in the review and except in the early stages of reaction, the 
as sumptio no f sample homogeneity is not likely to cause any 
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significant errors, 

(9) SEOOITDAEY BEAOTIOIIS ; The low feasibility of C-CO 2 
and CO-Og reactions has already been disci^.ssed in the thermo- 
dynamic review. Even if secondary changes do take place, a 
first ord.er rate expression wath respect to the concentration 
of oxygen is satisfactory enough for the empirical description 
of the gross surface changes, 

(10) SUBLIilATIOlT OF GRAPHITE .- At the low operating 
temperatures under consideration, the vapor pressLire of graphite 
is so negligible as to render the possibility of sublimation 
remote. 

With these assumptions in mind, evaluation of various 
models for determining the first order rate constant' from the 
experimental weight loss curves may be considered. The rate 
constants thus determined will be used in the main stream of 
analysis of the problem of a noncatalytic reaction in a pore. 

As mentioned in the experimental procedure, a system of a, stag- 
nant air column over the surface of a reacting graphite pellet 
is chosen for the determination of rate constant. Mathematically 
this would mean solving of one dimensional unsteady state diffu- 
sion equation wdth reacting surface boundary condition. The 
diffusion may be considered (a) in the gas phase (b) in the 
solid phase or (c) in both the phases, (c) would involve the 
solution of two simultaneous partial differential equations 
v/ith common interface boundary condition, (b) requires the 
kno¥;ledge of the effective diffusivity through the reacting 
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soa-id e„nd tills is not always availatile. Therefore physically 
and mathematically (a) becomes a convenient choice. Referring 
to the figure 4-2, the necessary diffusion equation and boundary 
conditions may be ivritten a,s 




AIR 

GOLUIM 


i i i i i 1 

OXYGEl 


ts^vwwws 


y=i 


y=o 


V Grapliite sample 
RIGUHE 4-2 

4-22 can be obtained as 


n 


R(s,y). = .) 
R(s,y) = 


k 


.3/2 


...4-19 


' a?' 

C(0,y) =1 

c(t,«> ) = 1 or C(t,l)=1 ..4-21 

D -||r (t,0) = kG(t,0) .. .4-22 

If a new variable f=1-C is intro- 
duced and further if the Laplace 
transform of f(t5y) is defined as 

CO 

(-yy) e“®%t. .. 

= RCsjy) ...4-25 

the solution in Laplace domain 

for the equation 4-19 with the 

boundary conditions 4-20, 4-21, 



f 71 s 

. :v 

k[R(s,0)-l/s] 


^ b^2J 


. . .4-24 


-fl+e~^^Ws 


.(y-SDx/Vf : 


J..4-25 


Note that the , solutions 4-24 and 4-2 5 are valid for infinite 
and finite length of the air column respectively. 

Inverse transforms of the terms on the right hand 
sides of the equations 4-25 and 4-24 Eire not readily available 
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in tlie standard tables and it is dcubtfiil whether the closed 
forms of the same are analytically possible. Efforts to 
obtain the haplace inverses and the general solution C(t,y), 
need not be persued, because in the determination of the 
surface roaetion rate constant, it is sufficient if the varia- 
tion of conce n t r ation with respect to time, only on the surface 
of the grarhite i.e, at y=Q, is known . So substituting y=0 
in the edimtioris 4-24 and 4-2'5, we get 


’(s,0) 


S( <.j 8 + lc/v5) 


, . , 4—26 


and 


E(s,0) = 


k(l 




1 


.4-27 


s [ VDs ( 1 +e-2 1 W s ) + 1 , ( ^ .e-2 ^^'^3 ) ] 


In the actual experiment, the length of the stagnant 
air column 1 is of the order of 8 cms, Sy interf eromotry [60a] 
it is well established that the significant gasphase concen- 
tration variations over the surface of the sample are known 
to occur within a distance of a few millimeters only, There- 
fore, for all practical purposes, the distance 8 ems is large 
enough to be considered as infinite, The inverse transform 
then need be considered only for the equation 4-26 so that 


f(t,0) = 1 - erfeCksft^) 

...4-28 

2 

c(t,0) = erfc (k N/t/B) 

. . .4-29 


Making use of the definition of the heterogeneous surface 
reaction rate constant, the rate of -loss of solid may be given as 
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= k erfc(kN/t7D) 


..4-30 


’wliere S is the reacting surface area and W is the Y;eight of 
the sample. Though in principle, the simultaneous determination 
of k and I) should be possible from the above equation, attempts 
to use the experimental weight loss curves for this purpose 
v/ere not fruitful in yielding consistent values of k and D, 
Therefore, evaluation of k from equation 4-30 was tried by 
using known values of D from literature, but again inconsistency 
of k persisted. The root cause of this inconsistency lies in 
the incompatibility between reliable experimental data and -' ■ 
accurate numerical values for error functions . In the Initial 
stages of reaction, i,e. for smaller va.lues of Icyt/Ii, the 
error function evaluation is accurate but the experimental weight 
loss curves are not reliable, and vice versa operates in the 
later stages of reaction. Even a special computer programme 
devised to estimate error functions accurate to the 28th decimal 
place (for the parameter k \jt/D ‘y A) did not prove to be of 


any use. 


let a 


another way to simplify the equation 4-30 for 


large values of kN.ft/1), may be to express the error function 
in the forra of’ an asymptotic expansion. Considering the first 
two terms in the series, the equation 4-30 may be written as 


1 ^ 
S dt 


C, VP 
Ao ^ 

sJtc t 


2k^t 


...4-31 


Integrating this equation, weight loss may be given as 
S C. \/F ^ ^ • 


W -W 
0 


sflT 


[ 2 ■/t^ + 


...4-32 


k V-t 
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Again, the weight loss cuttos could not fit in this equation 
for any consistent value of k. 

from the experimontal weight loss curves it can easily 
he seen that the time for which the curve is nonlinear is much 
less compared to the time for which it is linear. This suggests 
that a quasisteady state analysis of the equation 4--19 ciay he 
used with advantage in the determination of k without any 
significant loss of accuracy. The error involved in the quasi- 
steady state assumption may he ostimscted with the help of Wen’s 
chart [ 3 ]. In the present case, the initial oxygen concentration 
gm mol/cc, the solid voidage e dO. 3, the stoichio- 
metric number "a” in the reaction C+Op = COp is ONE, and the 
solids concentration gm, mol/cc. Therefore, the 

characteristic parameter (eOj^)/(a x lO"*^ and the 

corresponding error is less thsn OITE percent, Eor a quasi- 
steady state, the diffusion equation 4-19 then becomes 

.2 d’e . ■ - 

2^ = 0 or = 0 ...4-33 

2.7 dy" 

The hoiindary conditions ?/ith reforenco to figure 4-2 are given as 


D ^ (0) = Wyo) 
cpi) = 

The solution is = kC^^y/(lk+D) + ECj^^/(lk+E) 
As before, (-l/S) || = kDC^^^/(lk+D) 
low, -dW = ^ S dl 

1 = (W^-W + S 1^)/ S) 


. ..4-34 
. ..4-35 

... 4t36 

. ..4-37 
...4-38 

...4-39 
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Eliciina-ting 1 between the egr-ations 1-37 and 4--39 and integrating, 

' -i- ..,1-40 




2S if DC , 
Ao 


SD C 


AO 


and 


dW 


ii _ o L 1 

= ■trx’* — + tT 

S^f DC, ^ ^ 
js Ao 


Ac 




^ D ^ 


. .. 1-11 


■The plot of ^ versn,s (7/^-W) should in principle be able to. 
giTe the valne of k, bnt the intercept contains a tern in 1^ 
vdiich seriously affects the accuracy of k. There are a nunibei’ 
of empirical methods by which the weight loss ci^rves can be 
linearised but the tna thematic al loop holes are too many. Assum- 
ption of a complicated rate expression for the surface process 
has attracted many investigators but the numerical difficulties 
associated with the nulticonstant system out?;eigh the possible 
advantages of linoarisation. Considering all these difficulties 
it was decided advisable to fall back upon the conventional 
method of determining the rate constant from the limiting 
linear weight loss curves. A number of investigators have found 
this approach satisfactory [2], 

Before proposing, testing, and investigating the merits 
and demerits of vesrious mathematical models for the main problem 
of reaction in a pore, it may be advisable tc make a precise 
state'O'ient of the problem xmder consideration. If is the 
inner radius of the graphite c^^linder, is the external 

radius of. the graphite cylinder, C^^^ is the concentration of 
oxygen in the gaspliase at the mouth of the pore or the initial 
concentration of oxygen both inside and outside the hollo'w 



10 2 


space of the graphite cylinder, is the concentration of 
oxygen at any time and axial distance, i,e. the yariable 
concentration of oxygen inside the pore, R is the variable 
radial distance, z is the variable distance along the axis 
of the cylinder, t is the tine of reaction, 1 is the •uniform 
temperature of reaction, H is half the length of the cylinder, 

D is the diffusivity of oxygen in the gas phase reaction mixture, 
k is the first order rats constant for the surface reaction 
of graphite Y>rith oxygon expressed in the units of length/unit 
time, and further, 'if the external surface of the cylinder is 
insulated against chemical reaction and only the internal pore 
surface is allowed to react, how will the parameters and 
vary with respect to the axial distance z and the reaction time t 
and how best can this data be used in the design of heterogeneous 
reactors? 


Since it is a problem of Poiseuille diffusion in the 
pore space and a heterogeneous reaction on the pore surface, 
one is naturally inclined to consider the general equation of 
diffusion and solve it for‘ the react’ing v;all surface condition, 
T'he relevant isothermal unsteady state diffusion equation in 
cylindrical co-ordinate system is 

D ^ 




(R 


A ^ A 

^) + D ^ A 




... 4-42 


€>t “ R dR btTST 

formulation of the boundary condition for the reaction on the 
wall surface and the analytical solution of the above equation 
being difficult, the following simplified models have been 
considered. 
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(1) STOICHIOMETRIC MODEL * Ueglecting radial diffusion, 
the unsteady state material balance equation for simultaneous 
diffusion and chemical reaction of oxygen may be given as 


?lz 


Ri 


. . . 4-43 


This equation can easily be solved provided R^^ is a constant 
(i.e. independent of t and 2 ), but during the course of 
reaction R^ does vary, say according to a functional relation 


Rjl_ = g(t, S;) . ...4-44 

The nature of this ’’g*' being unknown, the equation 4-43 cannot 
be solved, and to add to the dilemma, "g" cannot be determined 
unless the solution to the equation 4-43 is available. There- 
fore, a straightforward material balance approach will not 
yield a solution even if numerical techniques are used, 


2. AVERAGE RADIUS MOREL ( CYLIH33RICAL PORE MODEL) ; In 


view of the paradigmatic value of the nature of variation of 
R^, it should be logical to make some assumption in this 
direction. If a mean radius R^ is assumed, the equation 4-43 


becomes 


D 


y Ca 


2k C , 
A 




• . . 4—45 


Defining a dimensionless axial distance x=s/H and a dimension- 
less concentration C = 


D C 2k C _ ^C 

R. 


...4-46 
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By tie method separation of variables the general solution 
of the above equation may be written as 

p r r— ■ 

C =-■ A^cosh{ jx '^/l- alnli(jxyi- B^) ...4-47 

i.,>' 

where g, A..J and B.^ are constants and the 

Modified Thiele Modulus 3 = lh^'^2k/(I)S^) ...4-48 

The relevant boundary conditions are 


c)C/ 


C(t, 0) = 1 

(t,. 1) = 0 
1 


dx' ’ 
C(0,x) 


4£(t, 0) = 2£:£ .o = 3^5 
BE, 


. ..4-49 

...4-50 

...4-51 

. ..4-52 


where 0 is the average concentration over an axial distance/ H 
corresponding to the average radius Rj 5 i.e. 


C = 




C dx 


pL 


=.rc 


dx 


...4-53 


% 

o 


Since it is difficult to evaluate the constants in 
equation 4-47 by accommodating four independent boundary 
conditions 4-49, 4-50, 4-51, and 4-52, the procedure for 
calculations will be indicated by assuming a solution of the 
type 


C = G(x, t, Ex, 0) 


...4-54 


Now, the change in average radius Rj_ can be related to S’ as 

Sr. 


■1 e 




4-55 


At the start of reaction, the distributionof concentration 
being uniform, C = 1 and therefore 
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v;liere (d,Rj_)-| is tlie ch^ge in average radius corresponding 
to the first time increment (A-t)^ 


k C 


/. (Rj_)^ = (R^)^ + (A\)^ = (R.)^ + ...4-57 

where (Rj_)Q and (Rj_)^ are the average radii at the beginning 
and after the first time interval respectively. 

Using this (Rj_)pin equation 4-46, the concentration 
distribution of the type 4-54 can be obtained, and by going 
through similar steps, (R^) 2 , the average radius after the 
time increment (/it) 2 j can be evaluated. The procedure can 
be repeated for further time increments till the pore mouth 
radius becomes equal to the external radius. However, once 
this stage is reached H and hence the modified Thiele modulus j 
assume different values for each incremental cycle. This 
•"forward march tectinique”, -with step variations in o, can be 
employed to determine the distribution of concentration and 
the average radius at various stages of reaction till the 
whole cylinder is almost completely reacted, , 

Obviously, there is a limitation to the utility of this 
method to the extent that at any time t during the course of 
reaction, it fails to give the distribution of the radius in 
the axial direction. This limitation can be overcome by 
formiilating separate sets of differential equations for a 
number of small grids each of an axial length Ax, and computing 
C and (Rj_)j^ separately for each of the grids. However, it may 
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be noted that if the total number of grids is 1, there will 
be an K fold increase in the number of computations. Consi- 
dering all these inadequacies and inaccuracies, it is doubtful 
whether the average radius model can provide a satisfactory 
mathematical description of a noncatalytic reaction in a 
pore. 


3. CONICAL POES MODEL S Having considered an average 




radius in, a cylindrical 
pore model, the next logi- 
cal step should be the 
assumption of an inclined 
pore wall. Figure 4-3 
illustrates a closed pore 
of a conical shape, the 
height of which is H and 


I’lGUES 4-3 


the pore mouth radius 

Let half the angle subtended 


at the apex be Q i.e. tan”^ j and for convenience, the 

origin is located at the centre. Keeping the rest of the notations 
as before, the material balance over the frustrum of an ele- 


mental cone of thickness £^..Z is given by 

'ho p p j 

(-D-^ Hh|) - -.0 




-kC ^ It ( R^+R^+ A R ) ( A z^+ A R? ) " 


+vt (R.+ AR.)2+y,tR2^(R.+ ARi)2V || 


...4-58 


Substit-uting ik tan9 and = z tan9 9 and taking the 
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limits as — 'y o , 




A 


c)z^ 


+ 


2 

z ^z 


2k 

irn 


( 1 + “ban ©) 




~W 


A 


...4-59 


for a quasi steady state, the above equation becomes 


where 


d 


dz' 


h' 


' z dz 

2 (2k H^) 

■ DR. 
im 



4 - 


tan^© ) 


. .. 4-60 

. .. 4-61 


The necessary boundary conditions are 

C^(H) = ...4-62 

^ c 

-3) (0) = kC^(O) ...4-63 

The solution of the equation 4-60 then becomes 


Ik (H)i h [ 

°Ao ^ I, [ 211 ] 

where represents the first order modified Bessel function 
of the first kind. Thus for a quasisteady state, the above eqn 
represents the concentration profile of oxygen along the length 
of the pore. It is also implied that the time required for 
the pore to change its dimensions significantly is large as 
compared to the time required to establish the quasisteady 
state concentration profile within the pore. The rate at which 
graphite reacts with oxygen may therefore be given by 


RATE 




[2hv4/?l 


[2h] 




4-65 
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After a finite time interval At, the pore radius increases 
hy A and the relevant material balance equation is 


A R. 




^ = RixlE or 




At 


, . , 4-66 


Equating the poremouth flnw per. unit surface area with the 
reaction rate within the pore, under a pseudo steady, state 
condition, 


RATE = 




iIo(2h) 


1 


H' 




... 4-67 


where is the zero order modified Bessel function of the 


first kind, Comhining the equations 4-66 and 4-67, 


At 


R. ^RC. h 

xtn Xi-O 


Po^2h) 


(l+tan^© )^ 1^1 


. . . 4-68 


It may be noted that the equation 4-68 makes the following 
implicit assumptions; 

(a) Small changes in radius do not significantly 
alter the apex angle of the cone, 

(b) The pore is long enough to keep the oxygen concen- 
tration at tho apex almost zero, 

(c) And as a corollary, the apex of, the cone does not 
shift significantly away from the origin. 

H, the parameter h is so large that the equation 
4-68 can be further simplified as 


At 


^Ao 

f H 

■'s 


(2kRj_D)'^ 


... 4-69 


In the initial stages of reaction, the pore geometry 
being cylindrical, equations 4-68 and 4-69 obtained from a 
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conical pore model, obviously cannot be applied. Moreover, 
even if the pore develops a conical shape during the course 
of reaction, equation 4-68 or 4-69 cannot accurately predict 
the wall profile because the analysis presented is that for a 
closed pore while in reality the pore is open at both ends i.e, 
a divergent-convergent-divergent pore. 

However, the equations 4-68 or 4-69, v/hich represent the 
rate of reaction inside a pore, may be used to formulate an 
index which gives the ratio of reaction rates on the internal 
and external surfaces of the solid. If the external surface 
of the cylinder were also to react simultaneously, the rate 


at which the external radius He would decrease may be given as 

k C, 

...4-70 


ARe 


^Ao 


At 




Combining equations 4-69 and 4-70, 
AE. ( 2kR.D)^ 

^ IL? 

ARe ” m " 


...4-71 


where 4> denotes a predict tv tty parameter". It can predict 
whether the burning of graphite takes place essentially on 
the internal or external surface depending on whether ^ > 1 
or < 1 respectively. It may be noted that in the equation 4-70 
the term can suitably be modified to accommodate the 
resistance due to convective mass transfer from the bulk of 
the gas phase to the solid surface. Standard Husselt number 
correlations for geometries other than a cylinder may be used 
to generalize the predictivity parameter. It is hoped that 
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this quantitative approach may prove to be useful in the 
optimum design of operating parameters in a heterogeneous 
reactor, ^ 

4. MITE OlfgEBENCE MODEL ; After exploring the possi- 
bilities of cylindrical and conical pore models, one could 
think in terms of various nonlinear geometries in their natural 
sequence. However, due to the mathematical complexity intro- 
duced by these, the analytical solutions of the equations 
become extremely difficult. Moreover, the assumption of one 
single nonlinear profile need not necessarily be valid for 
all tim-e intervals, temperatures and pressures of a single 
system - let alone different systems - Under the circumstances, 
the best way would be to consider a suitable numerical method 
to solve the simultaneous partial differential equations of 
the system. Irrespective of whether it is a forward march 
technique or an iterative technique, a stable finite difference 
scheme forms an essential ingredient of any numerical effort. 
Erom time to time, a number of new stable finite difference 
forms for quasilinear parabolic partial differential equations 
are reported in the literature, with a view to increase the 
size of the increments, and hence reduce the computation time 
required for the numerical solutions [ 108, 109, 110, 111, 112], 
It is probably not possible to make any absolute statement 
as to the '’best" finite difference scheme for solving a parti- 
cular partial differential equation. However, a good rule 
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appears -fco be. the following; given a maximum permissible 
error^' use that method which will give an approximate solutioi 
?;ithin the specified error with the least amount of machine 
time. It appears unlikely that any single method will be "best” 
from this point of viev/ even if the choice is restricted to 
parabolic equations only. It, therefore, becomes necessary 
to investigate a number of different finite difference schemes 
so that necessa.ry information is available when needed for the 
solution of a particular equation. 1 stability analysis with 
reference to lees* finite difference scheme [108] for the 
most general quasilinear parabolic differential equation, is 
illustrated belov;. It could as well bo modified and applied 
to other finite difference forms under consideration. 


(d) Milton Lees* Extrapolated Difference Scheme ; Let 
w(x,t) be a smooth solution of the quasilinear parabolic partial 
differential equation 


K(?c,t,w, -^) ...4-72 


3x 


in the rectangular region 0={(x,t)i Ol.x‘1.1, 0 <.t<.T} 


such that 

w(x,o) = (j)(x) , 0 < X < 1 ...4--73 

w(o,t) = $(x), 0 < t < T ,,,4-74 

w(l,t) = rfx), o < t < I ,..4-75 

where 4” >4"*" are specified functions. It is assumed that 

the functions K(x, t,p,q) and f(x,t,p,q) are defined and smooth 

2 2 

for (x,t)eQ and p +q < and that there is a positive 
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constant u sucli that, in this region, 

K(x,t,'p,q)> y >0 ...4-“76 

It follows from 4-“76 that the differential equation 4-72 is 
uniformly of parabolic type and, consequently, w(x, t) is the 
only solution of the boundary value problem given by the 
equations 4-72 to 4 - 75 . Incidentally, the equation 4-72 becomes 
"semilinear" if the functions K and f are independent of q. 

If h is a positive -integer and h >o is defined by the 
relation (l+l)h=1, then G(h) = {ihj i = 0, 1 ,... .11+1 },. . 4-77 
If 0 < : 21 g <1 then v/e let J be the largest integer such that 
kJ < 1, let C (h) be the real linear space of all functions 
u(x) defined on the finite set G(h). Obviously, ^(h) has 
dimension N+2, Ihe concern no\v is with the construction of a 
difference scheme for the boundary value problem 4-72 which, 
for all sufficiently small h and k, determines uniquely a 
sequence u^, u^ , U2,..,.Uj of functions from 5(h) such that 

(I) u = 4 on G(h) and, for each 3^1, the N numbeig u.(h), 

- J 

u^(2h), u.(Nh) are determined directly as the solution 

ti 

of a linear system of equations. 

(II) There is a 6> 0 and a constant 1> 0 , independent of h 

2 2 

and k, such that max ju .(x)-w(x, jk) } < A(h +k ) whenever 

xeff(h> 

/ {h^ + k^T^s. 

It is implied in (I) that the numbers u.(h), u.( 2 h), .... 

J ti 

,. .,U4(llh) need not be determined by an iterative process 
J 
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applied to nonlinear system of equations. Condition (II) 
states that the sequence u , u^ , furnishes a second 

order approximation to the solution of -the boundary value 
problem 4-72. 

It is Imown [ 110 ] that the standard Crank-Ficolson 
difference scheme for 4-72 satisfies condition (II). It is 
also known that it fails to satisfy condition (I), except in 
the case where 4-72 is a linear differential equation, which 
in the present context is out of consideration. In 1963, 

Douglas and Jones [ill] formulated two modifications of the 
standard Crank-lTicolson difference scheme for 4-72, which 
satisfy condition (I), IWien 4-72 is a semilinear equation, 
they proved that their predictor-corrector Crank-llicolson 
difference schemes also satisfy condition (II). In a general 
qu^silinear case, however, they obtained only the weaker error 
estimate 

]u.(±) - w(x,3k)j < A(h^ + ...4-78 

In 1967, Milton lees [1O8] put forward the forthcoming extra- 
polated scheme which, in a.ddition to satisfying the conditions (I) 
and (II) in the general quasilinear case, requires approxi- 
mately half the number of evaluations of the functions It and f . 
as the predictor corrector methods in [III]. However, it was 
found necessary to add to condition (II) a mild restriction 
is: :2Xh for some fixed positive constant^ , In 1970, H.B. 

Keller [109] has presented a new "box scheme" which is 
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’unconditionally stalle, simple and easy to programme, and 
efficientV However, this box scheme is excluded from the 
present consideration because it contains some newly defined 
derivatives which are useful only when the coefficient of 
2nd order derivative is a function of x and t. 

Now, coming to the extrapolated scheme, two completely 
inessential simplifications are made for notational convenience, 
the first one being the independence of K and f with respect 
to X and t, and the second one being the identical vanishing 
of and Thus it is supposed that w(x,t) satisfies 



3w 

3t 



+ 


f(w. 


Sw '' 

3x^' 


w(x,o) = (x) o<x< 1 

w(o,t) = w(l,t) =0 o < t < T 
Defining the standard difference operators 
the usual way, 



...4-79 
...4-80 
.. ..4-81 

1)_ and in 


h D^u(x) = u(x+h) - ...4-82 

h D_u(x) = u(x) - u(x-h) ...4-83 

2h DqU(x) = u(x+h)-"^(x-h) ...4-84 

Let C 5 (h) be the N-dimensional subspace of ?(h) consisting 
of those functions u ei; (h) for which u(o) = (uCl) = o. Let 


G(h) denote the intersection of 0(h) ?/ith the open interval 

p 2 i* 

0 <x <1, and let t = (h +k ). 

Two sequences Uq,u^ ...... ,u^ a^id u^,u^, , , . . ,Uj_i will 

now be constructed in ^^(h) by process of induction. Let 
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Now, assume tlat , , . . . ,u^ and uti , < » . » . ^'u ^ , i. have 

already been defined. Then u.-,- is defined to be the "unique 

D"T" I 

sol"ution of the linearized Crank-Nicolson difference e.quation 




+u.)+kf(u._^,D "u. i) .. 4“86 
3 3“^2 ^ J "*“2 


for X e ^(h), and u^_^^y '2 i-s defined directly by the linear 


extrapolation formula 

. _ 1 


j+3/2 - 2 ^' 3+1 " 2^3 


, . .4-87 


It is clear that the t'wo sequences <u .> and <u.,i> are well 

' • 3 2 +' S , 

defined, and the sequence <u •> satisfies condition (I). It 

J 

may be noted that the justification for the title "Extrapolated 

X^ifference Scheme" stems from the equation 4-87. Incidentally, 

if f is a linear function of p, .it is permissible to replace 

f(u^^^, u^_^i)by, a linear function 

of q, it is permissible to replace f(u-,i., 1 ' u..^) by :■ 

3 + 2 0 J-r2 

By making use of some properties of the linear space . 
^^(h) and by means of a straightforward but lengthy computation, 
based on Taylor's formula, Lees has proved the following 
relations for 


w_i = Uj^ + 


.. . 4-88 






”3+3/2 = I ”3+1 - + Ej+t 

“o”3+3/2 = K”3+1‘*®o'’'3+®J+* 


...4-89 


... 4-90 


... 4-91 
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= D Ui + R? 

O 2 Op- p- 


...4-92 


where, for 0 <3 < J, 

}IVi i 


<1 |E. il 1 +1 Ie'? i| ! < A. t‘ 

J+£ . <» 0+^ “ I 


...4-93 


Here is a constant independent of hjk:, and 3 , Without loss 
of generality, it may be assumed that 


i 1 


wj I < l^k 
^ 00 


. ..4-94 


If error functions = ^'/--Uj; and z-'i = w . , 1- u.^i are 

3 3 3 3+s . 3+'g' 3+2 

introduced, z. and belong to 5 (h) . It follows from 

J 0 '^ 2 O 

equations 4-85 to 4-94 that 
zi = Ri, = Ri 

S "2 0 

^3+3/2 = i ^3+1 ■ ^^3 ^ ^3+3/2 
^'o^ 3+3/2 = i^o^3+1 - i®o^3 ^ ^j+3/2 


. ..4-95 
,.. 4-96 
...4-97 


^‘^"^3+i’ V3+i^^"3+r"3'' =1 " 

^^~l^"3+it ^’3+l""^3-^"3+r^3^^^'^3+1 . 

^3+if"3+1-V^3+1-^3^? 

where the partial derivatives appearing in the above equations 
are evaluated at points called for by the mean value theorem of 


the calculus. 

for the solution v; of the .equation 4-79 > le'fc 

. m = “‘rEiwU : If I] ...4-99 

and let A2 denote the maximum of I'^l > l“fqi ’ l”lfN l“|f| 
over the closed bounded region determined by Ip] +iqU 2m. 

If k < Ah, then there is a-” 6> 0 such that t < 6 implies that 
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2A^ ax <1, 

2 crx <ni^ 

2A^ T<]i , 

2(A+A^ )t < 1 , 
where a = 5 + 4 and 


. ..4-100 
...4-101 
.. .4-102 
.. .4-103 


A 


21( y2 16{A +A2)V li 

v^- L + +1) J e ...4-104 


^i+i ' L ' f L < CTT 

^3^ 0 :S 3 < 1 


^ow^the conditions (I) and (II) will be satisfied if iiy induction 
i*t is pi?ov6d for x < £ ^ 

.,*4-105 
• . •4—106 

a Tlew of 4-93, 4-95, and 4-100 to 4-103 It follows that the 
inequalities In 4-105 and 4-106 hold for 3=0. Assuming that 
4 105 holds for all Integers less than or equal to 3, it will 
be proved that 4-105 also holdswlth j replaced 'by j+i . The 
inductive assumption and 4-101 Imply that the modulus of each 
partial derivative In 4-98 Is net_greater than Ag! Therefore, 

forming the Ig Inner product of 4.98 with 2(zj+i-Zj) and using 
4-76, v\?e get 

^"'^'"311 - "3- (-If) ^o"3tt) 

p 5 k 

j . - " ‘ 

O' aq 


“Z 




3+i '^"i+1-''3-("3i-1-^d>h + 2' (^3+1-",. it Wj^.-W. 


-'^3+1 - ^3^^ > 


...4-107 
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Prom 4-107, the following estimate may he given as 

+ [i l^j+4 I 1+ I I 1 z 1 it I |Zj+j ! 1+ I 1 1] 


+ i i Z;5+1“2^ M 1 i i ...4-108 

Prom the inductive assumption and from equation of 4-108 and 
from the general inequality relation 

2ah <:i, ea^ + - h^, e > q, • ...4-109 

another edtimate may he obtained as follows; ' 

[ 2 li -2A^ c T -2 e ] I I II + i ! z^_j_^l I ^ ^ I ! z^ | I ^ 

+k (A^4-A2) ^ C I I ^ ^ ^ ^ ^0^3+^ I ! ] + e "'l j I 

... 4-110 


how, choose 2e=u . Then, in view of 4-r102 and the inequality 


R 


3+1 


< A^ T 


it follows that 


llh+1 lb <il h II T, ^ V3+il0' 


+ 




D 


...4-112 


where A^ = 2 p (A^ + ^ 2 )“^ 

Since z = 0, the inequality 4-111 may he iterated to get 

° . 2 
3 2 2 _A. 


D - n=0 

,.,2 _2 


< (kA^A^ <3 + 2 y A^T) t + k -^■^[11 1 ^ 1 ^ 

... 4-113 ' 

where the relation 4-93 and 3k < Jk < T are made use of. ^ 

I I I I + 1 1 h I I < I 1 1 ^nl 1 -+ 1 1 V1I I + il Po^n 1 1 + 
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® view of 4-96 and 4-97, it follov;s 


from 4-113 that ^ A, + 1 6A,k E llz„l|^, where 

^ n=0 » ^ 


'D+1 

2 


D 


A^, = kA^A^ + 2;. u ^ A^ T + 4A^ A^ T ...4-114 


Using the discrete version of Gronwall’s lemma given in ref. [110], 

...4-115 


2 4 

z ... I , < A. T ' 

0+1 ]3 4 


el6A3l 


'D+1 


- /r 2 ^8A,f 

j+lllu 1 "4 ^ ® ’ 


It follows from 4-100, 4-102, 4-112 and 4-114 that 


(AJ^ < 


1 


+A. 


4- - /y 
which implies that 


[c^ + 21 ( + 1)] 


JL 
. 2 


. ..4-116 


...4-117 


< At; therefore 4-106 holds, 
with j replaced by 0+I . It follows from 4-96 and 4-97 that 


^0+1 


11^3+3/2!! ^1 |^o^0+3/2| I < 2 A_t^ + 2A^ T^+ 

+-|-| jU^z^l I < 2(A+A^) t + 4h“^ ...4-118 

^ince ^ I'l^i 11 ? hence by 4-103, we have 


z 


3+3/21! /■! 1^0^3+3/21 1 < 2(A+A^) + 4(1+ x^)h 


<t+4{1 +a^) t 


= (5+4X )t=0t »• .4-1 1 9 

which verifies 4-105, with 3 replaced by 3+I • Ihns the extra- 
polated Crank-Uicolson difference scheme given bj?- equations 4-85: 
4-86, and 4-87 satisfies the conditions (I) and (II). 

After ensuring the stability of the difference scheme, 

/ 2 2 

and after estimating an error bound of the order of V h +k , 
the basic differential equations 4-46 and 4-55 which describe 
the system of reaction in a cylindrical pore, may be recalled 



. ..4-120 


in the follo?/ing forms 


30 

3 1 



3^0 



2kG 

R 


D 

S 


as 

3 1 


kO^ = kC^C 


For numerical convenience 


...4-121 

introduce a new variable f=1-C so that 


- P 

St-= g2 



+ 


2k(l-f) 

___ 


“fi = ^ ^Ao 


. .. 4-122 


...4-123 


Writing these equations in terms of Lees’ extrapolated difference 
scheme given by equations ■ 4-85, 4-86 and 4-87, and denoting 
f(x,t) = f(n.m) and R(x. t) = R(n. m) . 


f(n,m+l) -f(n,ffl) = {f(n+1,m+l) 

2H^ 


2f (n,ni+1 )+f (n-1 jffl+l ) 


+f(n+1,m) - 2f(n5m) + f(n-1,tn)} 4-2k At(l-3/2f (n,m)+i-f (n,m-l) )/ 

(!■ R(n,m) -iR(n,m-l)) ...4-124 

kC, At , 

and R(n,m+l) - R(n,tn) = — r [l-^ f(n,m)+-if(n,m-1 )]. .4-125 

^s 

Standard computer programmes [l15]are available to solve the 
equations 4-124 and 4-125 simultaneously. The rank of the 
matrix of coefficients in the equations of the type 4-124, depends 
on the number of axial grids used in the computations. The 
method of calculations for ax = 0,2 i.e. for five axial grids 
over half the axial length H, is given belOTv'S 

The equation 4-124 contains terms at time level m+1, m 
and m-1 . To continue with the "forward march" the terms at time 
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level m+1 are expressed in terms of those at time level m and 
m-1. Referring to figure 4-4, the axial grid position n=0 



r- - P* 





© 

i 

1 .2 ,3 .4 ,5 

-6 .7 

4 

t i ! u 

1 1 


GRID POSIIIOHS OVER HAT. P THE 
PORE 

PIGUEB 4-4 

symmetry, the point n=7 forms the 
If (D At)/(H^ Ax^) = h 

(k it)/ Pg = a 

the equations 4-124 and 4-125 may 


represents the air space out- 
side the cylinder. The positions 
n=1 , 2, 3, 4, 5, 6 represent the five 
grids over half the length of 
the cylinder H. By virtue of 
mirror image of the point n=5. 

,..4-126 

...4-127 

he re-arranged as 


^ f (n+1 , m+1 ) -(h+1 ) f(n,m+l) + ^f(n-1,m+l) = -^f (n+ 1 , m) + (h-1 )f (n,m) 

f(n-1,m)-2k At(l- f(n,m)+i f (n,m- 1 ) )/(|R(n,m)-iR(n,m-l) ) 

... 4-128 

and B.(n,m+1) = R(n,m) + a (l- f(n,m)+i- f(n,m-l)) ...4-129 

If at time t=0, m is set at zero, we have 


f(n,-l) = f(n,0) ...4-170 

and R(n,-1) = R(n,0) ...4-131 

and f(n,l) and R(n,l) would then denote the properties at point 
n and after a time interval At. Thus it may be noted that if 
equation 4-128 is to be used to calculate the properties at time 


m+1 , in terms of those at time m and m— 1 , there would oe three 
unknowns at positions n+1, n, and n-1 . Referring to figure 4-4, 
the equation 4-128 may be written at positions 1,2, 3, 4, 5 and 6 to 
obtain six simultaneous linear algebraic equations containing 
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EIG-HT unknowns. However, in 'view of tlie fact that f(o,n 3 ) is 
always zero and f(5,m) is always equal to f(7,m) the number 
of unknowns reduces from EIGHT to SIX. The solution hence can 
be obtained by means of the following simultaneous equations, 

^ f(2,m+l) - (b+l)f(l ,m+l)+ | f(0,m+l) = - | f(2,m) + (b-1 )f(1,m) 

- ^ f(0,m) - (2k At) (1- f('!,m) + ^f {l ,m-1 ))/(^ R(1 ,m)-'|-R( 1 ,m-1 ) ) 

...4-132 

|f(3,m+l)-(b+l)f(2,m+l)+ ^f(l,m+l) = - Jf(3,m)+ (b-l)f(2,m)- 

f(l,m) -(2k At)(l- |f(2,m)+if(2,m-l))/(| R(2,m)- iR(2,m-l)) 

...4-133 

....4-134, 4-135, 4->-136, .and 

|f(7,m+l)-(lH-l)f(6,m+1) + |f(5,m+l) = - |f (7,m)+(b-1 )f (6,m) - 

^f(5,Gi)-(2k At)(l- I f(6,m)+ if(6,m-l))/(|-R(6,m)-iR(6,m-l)) 

...4-137 

There are no difficulties in the calculation of radii at various 
positions with the help of equation 4-129. Again, referring to 
figure 4-4, six independent equations for radii may be written 


as follows! 

R( 1,01+1) = R(l,m) +a (l-|-f(l,m) +if(l,m~l)) ... 4-1-38 

R(2,m+1) = R(2,m) + a (1- ^ f(2,m) + if(2,m-1)) ...4-139 

.... .4-140, 4-141, 4-142.. .and 

R(6,m-#i) = R(6,m) + a (l - |■f(6,m) + if(6,m-1)) ...4-143 


Obviously Lees* extrapolated difference scheme involves 
an implicit method of solving a set of simultaneous LINEAR 
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e«at,on. to Obtain the oonoentnatlons after eaob Internal of 
This IS a significant adrantage over other difference 
soheaes which Involve a set of nonlinear slanltaneous equations 
Which tahe greater aachlne tlae for nuaerleal solutions. 

Baae ^fference Sche;^ ; Dnser the constraints and 
=on itions mentioned In lees' extrapolated difference scheme, 

OBe [110] approximated the parabolic differential equation 

3x ...4-144 

-ans Of a famll, of following Implicit difference equations^ 
+ Uj + (1-S)I)^I)_U._^ = S(x,t,u I (D^u + 

■ -Ml"; ^ . 

^ analysis argument (to obtain =n i 

an I2 error estimate) 

similar to the one mentioned before p i. 

Delore, Rose showed that 

Nu ~ wj ( = 0 (W) 

qa ^ . ..4-146 

c Sin y ‘V’Q.l'u© of ^nid. + *\ i a 2 

gria latio x = k/h'^ provided that 

0 < 2(l-s)x < 1 /a 

where 0 < a,,< 


0 < a^< < a* 

*- aq - 


. ..4-146 


... 4—1 4 7 
...4-148 


Extrapolating from the known results of Tin 

linear equations with 

constant coefficients, Rose also con;jeotured that | |u-w | j ^ = Olk+h^) 
_ 7^ provided that 2s-l > o. It may be 

e that the conditions s=l and s=i reduce the relation 4-145 

° the standard l,^ CrsnhS^ difference schemes 
respectively. 
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He?/, coming to the problem of reaction in a single pore, 
equations 4-122 and 4-123 may be written in terms of Rosei’s 
finite difference scheme (4-145) as follows! 

f(n,m)-f(n,m-l) = s K (f(n+1,m) - 2f(n,m) + f(n-1,m)) 

H 

+ (n+1 , m-1 )-2f (n,m-1 ) + f(n-1,m-l)) + 2k At(l-f (n,m) )/ 

H Ax 

R(n,m) .,.4-149 

k At 

and R(n,.::iii) - R(n,m-l) = ^ (1-f(n5m)) ...4-150 

^s 

Eliminating the nonlinear source term E(n,m) between equations 
4-149 and 4-150 and using the same notations for a and b as in 
equations 4-126 and 4-127, 

-abs [f(n+1 ,m) f (n,m)] + (a+2abs)[f (n,ffi)f (n,m)]-abs [f (n-1 ,m)f(n,m)] 

+ (abs + bs]R(n,m-1 ) ){f (n+1 ,m)}+ , ( (-T-2bs)E(n,m-1 ) + (-ab+abs) 
f(n+1 ,m-l)+(-a+2ab-2abs) f(n,m-l) + (abs-ab)f(n-1,m-l) - a-2abs- 
2k At){f(n,m)} + (abs +■ bsR(n, m-1 ) ){f (n-1 ,m)}+ (a+R(n,m-1 ) ) (b(l-s) 
f(n+1,m-l) + (l-2b + 2bs) f(n,m-l) + b(l-s) f(n-1,m-l)) + 2k At 
= 0 ...4-151 

Note that the quantities ?/ithin the square brackets "[]" represent 
the nonlinear concentration terms and those within the curved 
brackets represent the linear concentration terms, both 

at time level m. The rest of the terras represent the properties 
at time level m-1. Thus the equation 4-151 has three unknowns 
at time level m, (one each at positions n+1 , n and n-1) which 
are to be solved in terms of known terms at time level m-1. The 
problem now is similar to the one represented by the equation 
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4-128, with a major difference of the presence of nonlinear terms 
in the eqaation 4-151. 

Again, referring to the figure 4-4, six simultaneous 
nonlinear equations of the type 4-151 may he written for axial 
positions n = 152 , 3 , 4,5 and 6 in the follo?;ing manner, 

-ahsi_f(2,m) f (l ,rn)] + (a+2ash)[ f (1 ,m) f (1 , n.)] -as'b[f (0,m) f(l,ni)] 
+(asb + hsR( 1 , m-1 ) ) {f (2 , m) }+ ((-1-2hs) R(l,m-1) + (ash-ah) 
f(2,m-l) + (2ah-a-2ash) f(l,ni-l) + (ash-ab) f (0, m-1 )-a-2asb-2k At) 
{f(l,m)}+ (asb + bsR(l,m-l)) {f(0,m))+ (a+R(l ,m-1 )) (b( 1-s)f (2,m-1 ) 


+( 1-2b+2bs)f ( 1 ,m-1 ) + b(l-s)f (0,m-1 ) ) + 21c At = 0 ,..4-152 

.1 ....... o .. o «. 4-1 53, 4-154, 4-155, 4-156.... and 


-asb[f (7, m)f ( 6, m)]+ (a.+2asb)[ f (6, m)f ( 6, m)j -asb[[f (5, m) f ( 6,m)] + 

(asb + bsE(6, m-1 ) ) {f (7,m) >+ ( (-1-2bs)R( 6, m-1 ) + (asb-ab)f (7,m-1 
{2ab-a-2asb) f(6,m-1) + (asb-ab) f ( 5,m-1 )-a-2asb-2kA t) { f ( 6,m) }+ 
(asb + bsR(6,m-l)){f(5,m) }+ (a+R(6,m-1 ) ) (b( 1 -s)f (7,in-1 ) + (l-2b+2bs) 
f(6,m-l)+b(1-s)f(5,ffl-l)) + 2k At = 0 ...4-157 

Again, noting the fact that at position n=0, f is always 
zero and the position at n=5 is the mirror image of that at n=7, 
the six unknowns at time level m can be evaluated from the above 
set of six nonlinear simultaneous algebraic equations, by means of 
a standard iterative technique given in reference [113], Once 
the concentration profile at time level m is known, the radial 
profile at time m can easily be evaluated by means of equation 
4-150, and the "forward march" with respect to reaction time can 
be .continued. Actual computations for s=0.3 and s=0.7 are 
presented in the chapter on calculations and conclusions. 
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(iii) Modified Backward Difference Scheme ; If w(x,t) 
is a four times boundedly differentiable solution of the para- 
bolic differential equation 


“li = J’(x,t,w, -||) ...4-158» 

in the rectangular region fl = { (x, t) ] 0 < z < 1 , 0 < t < 1} and if 
the notations assumed in lees' difference scheme are maintained, 
the approximation to wCx^t) is given by the solution of the 
difference equation 

= ^(x,t,u^_^, — ^)... 4-159 

where G. is the difference approximation of the smooth function 

J 

g(x,t) defined in the region , Following a similar line of 
argument as before, it may be shoTO that 


1 lu - w 1 L= 0 (k+h^) , ...4-160 

The difference equation 4-159 differs from the standard backward 

difference analogue of 4-158 in that the difference derivatives 

D, and ■|■(D,+D ) are evaluated u. rather than at u^. This 
+ + - 0~ I 3 

modified backward difference scheme has special significance in 
the present context because the parabolic equation 4-122 describing 
the reaction in a pore belongs to a class of functions F which 
can be expressed as 

F(x,t,w,p,q) = Ej(x,t,w,p) + F 2 (x,t,w,p)q ...4-161 

It is readily verified that in such cases the solution of the 
modified backward difference equation requires the inversion of 
linear system of algebraic equations. 



If the difference scheme 4 -irq o 

59 iw applied to the parti- 
cular problem of reacting -iv, 

equations 4-122 sni 

4-14:0 become 

+■ 1 ^ 

and E{n,m) = E(n,m- 1 ) + a( 1-f (n, m- 1 ) ) ’ ' 

Where a 4 b^oonstants gl.on by equations 4-126 and 4-127."' ^ 

The problem now is a si -n-ni -i -p,- ^ 

^ „ . sifflplified Tension of equation 4-128. 

following a procedure sitliar to „„ 

in Lees' difference 

SCilOtll©^ ’t/ll© C OTi P PD +• -i o 

entration and radial prof iles -mau Uo ■ 

r-i-uxij.es may be given as 

bf(2,„)-(2b+1)f(,^„) + ^ 21c At(i:f(l,a-1))/E(i,„.1) 

+f( 1,01-1) 

+f(2,m-l) 

... 4 - 1 6 ^ 

"““•"“'“•••-•--•••••. 4 - 166 , 4 - 167 , 4 - 168 . 
-(7,.)-(2bti)f(6,nHbf(,,.) . 2, ^ 

f(6,a-i) 

= ®'4.-1)ta(l-f(2,.-i)) 


0000*0* e oo*. a* 

" • • O « o o 


...4-172, 4-173,4-174. 


* o « O • • • * 0X1(3, 


B(6,ti) = H(6.m-1) + aCl-f(6,„-1) .. .4-175 

Ip comparison with Lees' difference scheme, computatlonrhere' 
are sampler, but the error hound Is higher. The oomputatlonal 
»ertts and demerits of the modified difference scheme will be 
aiscussed in the chapter on calculations and conclusions. 



127 


( iv) Explicit Forward Difference Scheme s Before going 
into tlie details of the scheme, it may be noted that no rigorous 
stability analysis is available for the explicit <iifference 
approximations of the nonlinear or seminlinear parabolic 
differential equations. However, on a number of occasions, [114-]? 
the explicit difference forms which are proved stable for the 
linear parabolic equations are extended to the nonlinear systems 
without any valid theoretical justification. If the actual 
numerical calculations are in reasonable agreement with those 
obtained from the standard stable implicit schemes, the extension 
is considered acceptable . The scheme is otherwise rejected. 
Moreover, no guide lines are available to decide whether a parti- 
cular difference scheme should or should not be extended to a 
particular kind of nonlinearity. This kind of approach, though 
not satisfactory, has been widely used in the solution of non- 
linear equations with considerable computational advantage. 

After having reconciled to this empirical procedure, scanning 
of various possible difference forms may be considered. Hov/ever, 
their number being too large and the numerical advantage of one 
modification over the other being uncertain and marginal, it 
may be advisable to start with the simplest of the explicit 
finite difference schemes. 

The basic differential equations 4-122 and 4-123 describing 
the reaction in a cylindrical pore may thus be written in the 
following difference forms; 
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f‘(n,m+l)-f (n,m) _ f (n+1., m)-2f (n, m) + f(n-1,m) 

At 


+ 


' Ax" 

2k(l-|-(f (n,ni+1 ) + f (n, m) )/(i(E(n, m+1 )+R(n, m) ) ) 

kC 


.4-176 


and ^ -^( i.i(f (n, m+1 )+f (n, >n) ) ) .. .4-177 


Eliminating the nonlinear source term R(n,m+-1) between equations 
4-176 and 4-177 and rearranging the terms as before, the foll07\f- 
ing quadratic equation can be written for f(n,Di+l). 


2 

af (n,ni+1) -[4R(n,m) + 2a+ab(f (n+1 , m) - 2f(n,m) + f(n-1,m)) + 

4k At] f(n,m+l) +{ f(n,m) + b(f(n+1,m) - 2f(n,m) + f(n-1,m))3- 
{4R(n,m) -i- 2a - af(n,m)} + 4k At{2-f(n,!n)} = 0 ...4-178 

Denoating the coefficient of f(n,m+l) by -b’ (n+1 ,n,n-1 ,m) and 
the constant terms by c ’ (n+1 ,n,n-1 jm) the above eauation may be re- 
written as 


2 

af (n,m,+l) ~ {b* (n+1 ,n,n-1 ,m) } f(n,m+1) + c ' (n+1 ,n,n-1 ,m) = 0. 

-...4-179 

To obtain the concentration at time level m+1 in terms of those 
at time level m, the equation 4-179 may be solved as 


f(n,m+l) = (b’ (n+1 ,n,n-1 ,m) +V\b’ (n+1 ,n,n-1 ,m) )^-4ac ’ (n+1 ,n,n-1 ,m) )/ 
(2a) ...4-180 

and R(n,m+l) = R(n,m) + ^ (2-f(n,m+l) - f(n,m)) ...4-181 

Referring to figure 4-4, the concentrations at positions n=1,2,3, 

4,5 and 6 can now be evaluated by six independent algebraic 
equations of the type 4-180 . This is precisely the advantage of 
an explicit scheme v/hich avoi(3s the solution of a set of simultaneous 
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liriear or nonlinear equations prevalent- in an impLlcit scheme . 
Once the concentration distrihution at time level m+1 is known, 
the radii at various positions can easily he calculated hy means 
of the equations of the type 4-181 and the" forward march" can he 
continued. 

Since the stability of the above difference scheme is 
dependent on the stability of the corresponding difference scheme 
for the linear parabolic system, it may be worth analysing the 
truncation and round-off errors and the constraints placed on 
the variation of "radius of convergence" in -the following linear 
parabolic differential eauation 


-3E = 1 
9^ 


w 


.. .4-182 


9x" 


Discretising the variables in terms of lay lor series,[ll4] 


-w(n,m+l) - w(n.m ) 9 w^ .. i 

At . “ 21 




o o « • « 


9t 


21 


...4-183 


w(n-fl.m) - 2w(n,m) + w(n-l.m) 

,2 ^ ^ 
ax 

If u(n,m) is the difference approximation of the exact solution 
v>?(n 5 m), the explicit forward difference scheme for the equation 
4-182 may be written as 

u(n,m+l) = u(n,m) + — (u(n+1,m) - 2u(n,m) + u(n-1 ,m) ).. 4-184 
To start with, von Meumann test may be applied to obtain the 
criterion of stability of the scheme 4-184. If i,9, and X 
respectively denote /-i', a real number and a complex number, the 
the exponential perturbation expression 
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u(n,m) = ...4-185 

may be substituted in 4-184 and each terra may be divided by 
einO giraX obtain X as tbe root of the equation 

s' = 1-4 ( — ~ 2 ^ sin^ (^) ...4—186 

Ax 

The stability criterion requires that the imaginary part of X 
be nonnegative, or in other words that for every real 9, 

‘ e-^* - 1 ...4-187 

This is the case if and only if 

( At)/( Ax^) < l/2 ,.,4-188 


(Note that the corresponding constraint for the nonoscillatory 
decay in the Crank Nicolson scheme wrould be . 


p —1 

( A't)/( AX ) ^{2(2(l-cos(s Ax) ) )> ,..4-188a) 

Under the constraint 4-188 concerning the ratio of the mesh sizes, 
not only should the round off errors generated by the finite 
difference scheme 4-184 have a stable behaviour, but also the 
approximate solution u(n,m) should converge towards the exact 
solution w(n,m) as ax-j-0. To prove the same, equation 4-184 may 
be rearranged as 


u(n,m+1 ) 


= u(n+1 ,m) 

A X 


+ ( 1 - 


-- ^0 u(n,m) + 
AX 


u(n-1 ,m) 

AX 

. ..4-189 


so that u(n,m+l) is expressed as an average of u(n+1,m), u(n,m) 
and u(n-1,m) with nonnegative weighting factors. According to 
Taylor's theorem ¥/ith remainder, the differential equation 4-182 
provides an analogous relation 
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w(n,m+l) = —^2 v/(n+1,m) -f (l- — w(n,m) + — ^w(n-1,m) 

+ 0( At Ax'^) ...4-190 

Subtracting 4-189 from 4-190, and introducing error term 

z(n,m) = w(n,m) - u(n,m) ...4-191> 

the difference equation for z(n,m) may be v/ritten as 

z(n,m+l) = z(n4-1,m) + (l- - ■ z(n,m) + z(n-1,m) 

Ax 2 Ax'^ hx. 

+ 0( At Lx^) ...4-192 

Noting the fact that for an initial value problem, 

z(n,0) =0 , . .4-193, 

the difference equation 4-192 together with the constraint 4-188 


gives 2 

Vi < % + ® ... 4-194 

where sup jz(n,m)[ ...4-195, 

^ n 

and G is some sufficiently large positive constant. Moreover, 
the condition 4-193 implies that z^ = 0, whence from 4-194, the 


basic inequality may be derived as 


m - 


m At Ax 


...4-196 


Also if rounding-off errors occur because only 1 decimal places 


are used in computing an actual numerical solution of the difference 
equation 4-189, then an additional term of the form 10 ought 


to be included on the right had side of the in equality 4-194. 
Therefore the final estimate of the total error becomes 

Zj^ < 0 m At ( Ax^ 4- 10“^ At"”*) ...4-197 

The result 4-196 implies that for any mesh ratio in the range 
4-188 the solution u(n,m) of the finite difference approximation 
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4-184 to the basic initial value problem 4-182 converges towards 
the exact answer v/(n,ni) in the limit as 4x -»■ 0, However, 4-197' 
indicates that in order to retain control over the round-off 
errors generated by a numerical calculation of u(n,m), choice 
of an unduly small time increment ( 4t) should be avoided. In 
this connection, it may be emphasized that the relationship 
between estimates of the truncation and of the round-off errors 
is closely interlinked despite the diversity of the sources from 
which they, arise. Thus in the ultimate analysis, the advantage 
of independent equations in an explicit form is countered by 
the restriction placed on the choice of mesh ratio and on the other 
hand, the advantage of unrestricted mesh ratio is an implicit 
form is countered b.y the occurance of simultaneous equations . 

Thus in a particular problem, it becomes necessary to investigate 
both explicit and implicit forms and then choose the best one. 

(v) Modified Richardson Difference Scheme ; Even though 
Richardson’s "unconditionally unstable" difference scheme which 
uses a central difference for the first order time derivative, 
is only of historic importance [ 114 ], a number of modifications 
over this classic scheme have been suggested from time to time[1l5] 
and some of them are stable too. In the present context, however, 
only one important modification proposed by DuPort and Erankel [11 6] 
is considered for investigation. 

Retaining the same notations as in eouations 4-182 and 
4-184, the original Richardson’s scheme for the parabolic 
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differential equation - 4-182 may "be written as 


u(n.[n+l) - u(n.m-1) _ u(n+1,m) - 2u(n,m) + ^(n-l.m) 4 -I 98 

2 A "fc ,2 • • • 

Ax 

Dufort and frankel [ 116 ] modified the above scheme by using the 
relation 2u(n,m) = u(n,m+1) + u(n,m-l) ...4-199> 

and obtained the following difference form for 4-182, 


u(n,(n+l) - u(n,m-l) _ u(n-i-1,m) - u(n,m+1 )-u(n,m-T)-ru(n-1 ,ai) 

2 At “ .2 • . or\n 

The advantage of this . modification is that it is an explicit 
scheme which is convergent and stable for all finite values of 
the mesh ratio At/ Ax^ . Forsythe and Wasow [ II 4 ] have given 
an elaborate stability analysis and have concluded that the discre- 
tization error involved in the approximation 4-200. is of the 
order 0( Ax ) for any fixed value of the mesh ratio At/ Ax , 
However, like in any other explicit form, the proof has not 
been completely; extended to the nonlinear cases. If the parabolic 
equation is of the more -general type (4-122 coupled with 4-123), 
at other places too, the termiu(n,m) should be replaced by the 
average ■ 2 '(i/n,m+ 1 ) + u(n,m-l))[ II4]. As in Eichardson's method 
the value of u after the first increment of time i.e. u(n,l) must 
be calculated by some independent procedure before the calculation 
by equation 4-200 can be started. 

Coming to the basic problem of simultaneous diffusion 
and chemical reaction in a pore, the modified Richardson scheme 
4-200 may be applied to the system of equations 4-122 and 4-123 
to obtain the following difference formas 
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f(n,m+l) - f(n,ni-l) = ^■5-^2 i£{n+^,m) - f(n5m+1) - f(n,m-l) + 

f(n-1,tn)) + 4k At(2-f (n,m+1 )-f (n,m-1 ) )/(R(n, tn+l )+R(n,m-1 ) ) 

... 4-201 

kC, At 

and R(n,ni+1) -R(n,m-1) = — (2-f(n,in+1 )-f(ri,tn-1 )).. .4-202 

Eliminating the nonlinear source term E(n,m+l) from the equations 
4-201 and 4-202 and retaining the same notations for a and h as 
before, the following quadratic equation may bo obtained for 
f (n,in+1 ) s 

(a+2ab) [f (n,m+1 ) ] ^ - {2ab f (n+1 ,m)-rfc2ab f(n-1,ni) + (2+4b) R(n,m-1) 
-4ab f(n,m-l) + 2a + 4ab + 4k At }f(n,Ki+l) + {2b f(n+1,ni) + 

2b f(n-1,m) + (l-2b) f(n,m-l)} {2R(n,m-l) + 2a-af (n,m-1 )} - 
4k At f(n,m-l) + 8k At = 0 ...4-203 


As before,- denoting the coefficient of fvnjm-tl) by -B' (n+1 ,n,n-1 , 
m,[n-l) and the constant terms by C * (n+1 ,n,n-1 ,m,m-1 ) , the 
equation 4-203 may be rev/ritten as 

a[f (n, m+1 ) ]^- (B' (n+1 ,n,n-1 ,m, m-1 ) ) f(n,m+l) + C ' (n+1 ,n,n-1 ,m,m-1 ) 

= 0 ^ ...4-204 

1 2 
f(n,m+l) = 2^ {B' (n+1 ,n,n-1 ,m,m-1) + ( (B’ (n+1 ,n,n-1 ,ni,m-1 ) ) 

-4a C' (n+1 ,n,n-1 ,m,m-1 )')^'> ,, .4-205 


Again, the concentrations at positions n=1,2,3,4,5, and 6 may be 
evaluated by six independent equations of the type 4-205. Once 
the concentrations at time level m+1 are known, the radii at 
various positions can easily be evaluated by means of equation 
4-202, The procedure is similar to the one used in the explicit 
forward difference scheme except that the term f(n, m) does not 
appear in the expression for f(n,ra+l). In fact, it is this 
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difference wliicli renders the scheme numerically stable. 

Presuming that the distribution of concentration and 
radii at any time is obtained by the use of one of the five finite 
difference schemes mentioned above, the problem then would be 
to formulate a model to determine the weight loss during the 


course of reaction. Por this 
the radii R(n,m) and ECn+ljm) 



PIGUEEt 4-5 


a linear wall profile between 
is assumed. Referring to the 
figure 4-5, the weight loss over 
an elemental cylinder of thickness 
Az is given by the weight of the 
elemental frustrum of a cone of 
thickness Az enclosed between the 
radii R(n,m) and R(n+ 1 ,m )5 minus 
the weight of the elemental cylinder 
of thickness Az and of constant;. 
radius R(n,o), the internal radius 
of the cylinder at time t = 0. 


P Az ' 

^ '^(R(n,m)) + R(n,m) R(n+1,m) + (R(n+1,m))^} 




-It (R(n,0))“ ^z 
s 

Summing these "AW"s over half the length of the pore H, 


.4-206 


ILW 



PgTt Az p p 

^ HR(n,m)) + R{n,m) R(n+1,m) + (R(n+1,m)) } 

5 ? 

^ PgTt AZ (E(n,0))"' ...4-207 

n=1 


Since the reaction is symmetrical about the centre of the pore, 
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the total loss over the entire pore is given by 

2 Az 5 , 2 

(W-W ) = 1 2 {(R(n,in))^ + E(n,m) R(n+1,Qi) +(R(n+1,ni)) } 

^ n=1 

-2 p 7iH(R(n,o))^ ...4-208 

s 

As will be seen in the conclusions, the equation 4-208 
gives a fairly accT;irate estimate of the weight losses at higher 
temperatures of reaction, i.e. when the wall profiles are linear. 

At low temperatures, however, a larger number of axial grids 
have to be considered to obtain more acciirate results from the 
equation 4-208. Graphical interpolation over the curved portion 
of the wall profile, serves as a convenient method of increasing 
the number of axial grids. 

■\j?hile investigating the finite difference model, it is 
important and interesting to attempt a mathematical model for a 
general pth order non-catalytic reaction in a pore, on the lines 
similar to those followed by Mehta and Aris [43] for catalytic 
reactions. The basic differential equations (4-120, 4-121) for 
a pth order reaction in a cylindrical pore would then take the 
following forms s 


ilG 

3t 


■R' 


^ 3^0 
2 




2k ’CW"”' 

Ao 

R 


. . .4-209 


Pg 


al ...4-210 

where k' is the pth order heterogeneous surface reaction rate 
constant having the dimensions of (moles)"^"^ ( length) time)”"' 
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A typical implicit scheme of lees applied to the equations 
4-209 and 4-210, gives the follovving relations? 

C(n,m+l) - C(n,m) = , {c(n+1,m+l) - 2C(n,m+l) + C(n-1,ni+l) + 

2r AX 

C(n+1,m) - 2C(n,m) + C(n-1,m)} -2k’ 0^^'' { | C(n,!n)-i C(n,m-1)} / 

R(n,ni) - tR(n,m-l)} ...4-211,and 

Ir* " 

a(n,m+l) - R(n,m) = { 2 . c(n,m)-i C(n,m-l)> ...4-212 

s 

The method of computations no?/ would he identical to the one 
followed for the equations 4-124 and 4-125, in Milton lees’ 
extrapolated difference scheme. The only difference is the 
presence of pth order terms in the constants of the simultaneous 
equations which alter the shape of the vmll profiles. 

On the other hand, if a typical explicit forward difference 
scheme is used, the equations 4-209 and 4-210 become 

C(n,m+1) - C(n,m) = {C(n+1,ni) ~2C(n,m) + C(n-1,m)} - 

.r Ax^ 

k' At(2)^ ^ {C(n,m+1) + C(n,m)} /{R(n,m+1) + R(n,m)} 

...4-213 

C T) 

and R(n,tn+1) - R(n,m) = {c ( n , m+1 ) + C(n,m)} .4-214 

^s 

As before, if an attempt is made to eliminate the nonlinear 
source term R(n,m+1) between 4-214 and 4-213, it results in a 
(p+l)th degree polynomial in C(n,m+1), Computations obviously 
would be complicated and inspite of the independent equations 
the explicit scheme may not have an edge over the implicit scheme. 
It would be of academic interest to note that a half order 
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reaction world give a 3/2 degree equation in C{n,ni+1) wliicli 
can be simplified to obtain a cubic equation in C(nsm+1). In 
the present context, however, no half order reaction between 
graphite and oxygen is contemplated and therefore the precise 
form of the cubic equation for various explicit difference 
schemes need not be worked out. It goes without saying that 
further calculations for the effectiveness factor for the pth 
order reaction would be still more complicated. 

The whole exercise in mathematical modelling may appear 
meaningless if reference is not made to the possible use of 
weight loss and profile data in heterogeneous systems. Attempts 
have been mads to iise single pore profile data for building a 
model to predict the overall rate of reaction on the inner and 
outer surfaces of a lump of a heterogeneous mass of ..any irregular 
shape and size. If Am is the maximum cross sectional area of the 
solid lump, and if Ro T r is the average radius of the internal pores , 
the given lump of mass ma.v be imagined to be a cylindrical pellet 
with right circular area Ap^ and with thickness 2H so that the 
volume of the imaginar.v pellet equals the original volume of the 
solid lump . The voidage is so adjusted that the volume of all 
the, pores (of radii Rqtt) imagined to be running parallel to the 
of the c.ylindrical pellet equals the original void volume 
of the solid lump . The overall rate of reaction of the solid 
lump may then be expressed in terms of the rate at v/hich the 
imaginary axial length of the pellet recedes during the course of 
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reaction. Bbweter, the ^ajor difficulty m such a scheme, is 
the absence of an analytical expression for the rate of reduction 
length of the pellet, tlumerloal data alone, may not bo 

complete validity of this new approach. 

nay be Interesting to obtain approximate solutions by method 

Of residuals and test the validity nf tvo 

vdiiaity of tne above parallel pore 

model. 


♦ • O e # 



5. CILCUTLATIONS A-ND COICLUSIOHS 

The accuracy of the calculations and the validity of 
the conclusions are primarily dependent upon the choice of 
realistic values of rate constant, diffusivity and density, 
and upon the judicious selection of the time increment At, 
the axial grid Ax and the radius of convergence b. Since 
the variation of the overall rate of reaction with temperature 
in a porous media or in a single pore is closely interlinked 
with the changes in rate constant, diffusivity and density, 
any error in these is' hound to have a significant effect on 
the quantitative predictions. Por instance, increased tempera- i 
ture has a positive effect on k and I) so as to enhance the 
overall rate of reaction while the same change in temperature 
has a tendency to decrease the density and hence the overall 
rate. Even the extent of relative variations in k and D have a 
profound effect on the establishment of unsteady state concen- 
tration gradients in the direction of axis of pore and hence 
of the axial diffnision fluxes. Therefore, as far as possible, 
attempts are made to select the experimentally determined values 
of density, diffusivity and rate constant. 

Reference [117] gives a fairly accurate and exhaustive 
data on the density of air and on the binaiy diffusivity of 
oxygen in air, at temperatures upto 5000°K. Theoretically 
predicted values are compared with the best available experimental 
data, and the errors reported are not more than 2,4^, Selected 
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values of densities at temperature between 873°K and 2073°K 
are indicated in 'Table 5-1. The contribution due to disso- 
ciation of oxjgen at temperatures above 1400°K has already 
been included. Tabulated values of diffusivity at various 
temperatures a:re not directly available in literature. However, 
reference [117] gives the thermal conductivity, specific heat 
and density and their dependence on temperature so that diffu- 
sivity can be calculated frora the ideal gas equation 


1 


Ic 


a 




...5-1 


The diffusivity values obtained from the equation 5-1 are 

slightly higher than the experimental values not so much because 

of the assumption of ideality but because the self diffusivity 

of air has been equated with the binary diffusivity of oxygen 

in air. One way of overcoming this flow would be to calciilate 

the binary diffusivity from the individual self-diffusivities, 

3/2 

b’ut a more convenient method would be to follow a ”7"^^ law" 
which is adequately accurate [106] v^ithin the range of tempera- 
ture under consideration. In fact this is equivalent to assuming 
an ideal -behavior at high temperatures ?/hich .is quite valid. 

The actual numerical values obtained from "T^^^ law" are given 
in Table 5-1, The base value of diffusivity at room temperature 
is taken from reference [105]. 

Coming to the computation of rate constant from the 
experimental weight loss curves, it has already been pointed 
out in the chapter on mathematical modelling that consistent 
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values of k Gannet be obtained from the fundamental equations 


of tbe type 4-30, 4-32, 4-41, and A-2 a Since the weight loss 
curves at temperatures 600°, 670°, 730° and 746°C are linear , 
^it is convenient to calculate the first order heterogeneous 
rate constant by the method of oonstait slope using the equation 




10 -^ 

2)(^q) 


= k C 


a « 4 5^2 


where S = reaction surface area Square centimeters, 

¥ = weight of the reacting Sample milligrammes, 
t ^ time of reaction seconds, 

4 initial concentration of oxygen in air St molesycm#^ 
The numerical values of rate ccaastants thus calculated are 


indicated on the respective graphs (page number 79 to 82)* 

At temperatures 870° and 1000° C the limiting slope method has 
been used on account of the presence of obvious nonlinearities 
in the initial portions of the weight loss curves. This method 
too uses the equation 5-2 in the calculations. 


Attempts were made to obtain k at high temperatures 
from the reliable experimental data at low temperatures i,e* 
experimentally determined values at 600°C, 670°C and 746°C 
were used to calculate the preexponential factor and activation 
energy which in turn were used to predict k’s at temperatures 
870°G and 1000°C, The predicted and experimental values differed 
hy a few hundred percent. This could be due to a number of 
reasons as discussed in the review of literature. Therefore, 
wherever possible, the rate constants at high temperatures 
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are taken from literature or are ca.Lcu.lated by judicious' 
extrapolation. Even though several numerical values of k 
(determined experimentally) are reported in literature, it 
is rather difficult to establish a precise scientific method 
for the selection of k for the particular variety of graphite 
used in the present investigation. Therefore the conclusions 
dra'wn hereafter should be viewed as applicable to a system 
with the properties montioned in Tabic 5-1 rather than to a 
system of graphite and air used in the experiment, 

Eesides the concept of a first older surface rate 
constant, there are other methods which enable a proper co- 
rrelation of the rate constant with temperature and porosity. 
The rate constfant (k^) may be expressed on the basis of 
reaction per unit volume of the porous solid, Levenspeil 
model [ 36 ] makes an effective use of k^ to describe a chemical 
reaction controlled regime, but iiie same does not hold true 
for a diffusion controlled regime. Moreover, the definition 
of k^ introduces extra terms describing the geometry of the 
solid in the basic equations 4-120 and 4-121, This increases 
the computation time considerably. At times, attempts are made 
to correlate k with porosity or percentage burn off of the 
solid. This, however, is an over simplification of the actual 
process of change in porosity during the course of reaction. 
Yet another way of incorporating the porous medium reactivity 
into the surface rate constant is to model the reaction in a 
macro pore by considering the axial', diffusion in the gas phase 
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and tlie radial diffiasion in the solid phase, This not only 
involres the determination of equivalent diffusivity in solid, 
hut also increases the computation time ~ even greater than 
the concept of k^. Since at high temperatures it is only 
diffusion vhich counts, attempts can he made to plot log k 
versus log T so that the slope of the curve which is approxi- 
mately equal to the index of variation of I) ^vith respect to 
temperaturE j can predict the relation between k and T. This 
again has its limitations in the unsteady state. And lastly, 
separate equations may be formulated for different ranges of 
temperature. This obviously has limited utility because of 
large number of constants involved in the description of the 
system, and extrapolation could he highly erroneous. Thus 
none of the above methods is really satisfactory for compu- 
tational adoption in the present system. 

It is of importance to note that, at high temperatures, 
there is no significant change in k because the reaction lies 
well in the diffusion controlled regime where k is adultrated 
by h. However, the present interest is unquestionably in an 
unadultrated k since the entire effort is to isolate the 

effect of k from that of D with special reference to the 

hyp o the t i 0 allj’’ 

reaction in a pore. It, therefore, become necessary to 
consider a nonporous solid at high temperatures and hence the 
exponentially increasing rate constant. To this extent, even 
at low temperatures, the porous graphite used in the experiment 
is equivalent to a nonporous solid with the surface rate constant 
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determined by the v/eight loss curves. The actual values of 
k used in the numerical calculations are indicated in Table 5-1. 
In fact, it would be ideal to undertake an elaborate investi- 
gation by considering all conceivable values of k for one 
particular value of I) and vice-versa. 

Coming to the model evaluation for the reaction in a 
pore, various finite difference schemes have to be tested 
with special reference to the stability, accuracy and the time 
of computation. In this connection, it may be noted tha,t 
even though tvro finite difference schemes for the two equations 
4-122 and 4-123 are independently proved to be stable, the 
stability of the simultaneous solution is not guaranteed [118], 
Moreover, there are no fixed rules as to which ofthe two 
equations should determine the choice of the grid. The boundary 
conditions influence this decision. In the present system, 
however, the nonlinear term being present in the parabolic 
equation, the whole stability analysis and the choice of the 
grid is centred around the equation 4-122, 

Besides the common problems associated with the finite 
difference approximation for a set of simultaneous partial 
differential equations, there is an additional uncertainty 
involved in the system of equations 4-122 and 4-123 regarding 
the time varying wall profile. Although equation 4-122 is 
formulated for a strictly cylindrical geometry, it is used 
during the entire course of reaction inspite of the fact that 
the geometry of the pore does not and cannot remain cylindrical 
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tlirougliout, Sucli usage is valid only when the profile of 
rile type given in figure 5-1 A can he approximated by the 
curved surfaces of a number of ri^t circular cylinders indi- 
cated in Figure 5-1B. Equation 4-122 then becomes applicable 
to each right circular cylinder separately. Ihe stability 



of the finite difference scheme thus becomes dependent upon 
the stability of the above -approximation which can be justified 
by the actual numerical calculations alone. 

Among the various finite difference schemes discussed, 
it is only in the explicit forward difference scheme that the 
size of the mesh ratio is constrained by the criterion of '■ 
stability. In the rest of the schemes, practically any mesh 
ratio is permissible provided the accuracy can be compromised. 
Thus it would be advisable to search for the guide lines to 
first decide the mesh ratio for the explicit forward difference 
scheme and then test the accuracy of other schemes by utilizing 
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tile same mesli ratio and mesli sizes. 

According to Keller [109], wiicn a nonlinear parabolic 
differential equation of the type 4-122 is coupled with a 
first order linear partial differential equation of the type 
4-123, it is preferable to choose the distance grid Ax 
first and then decide the value of time increment At accord- 
ing to the rules of stability. The axial grid Ax for an 
explicit forward difference scheme is chosen on the basis of 
the following considerations, 

(i) Jon liuinann's Stability Criterion t The constraint 
4-188 applied to equation 4-122 gives the condition 

b = (D At/XH^ Ax^) i.1/2 ...5-3 

Since At and Ax are two unknowns in 5-3, one more relation 
for Ax is necessary. Uormally, thco boundary conditions form 
the qualitative guide lines for the choice of Ax, However, 
no quantitative relations are available for the particular 
system under cons iders.t ton. Therefore, as a starting point, 
the condition for nonoscillatory decay (5-4) for the Crank 
riicolson scheme mey be borrowed -and applied to evaluate Ax 
in the explicit forward difference schorae. Therefore, from 

4-1S8a, ^ _ 

^ < [gin( 1C ...5-4 

li Ax . - . : 

Assuming a critical value of radius of convergence (b=l/2), 
the equations 5">f3 and 5-4 can be used to determine the low?er 
limit for Ax as follows; 

Ax >1/6 


. . . 5-5 
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The restriction 5-5 is normally applicable to a system in 
which partial derivatives of f do not occur in the boundary 
conditions s but if they do occur the number of unknovvns at 
tiniG level m+1 in equation 4-176 should not exceed the number 
of axial grids. There being only one unknowa is 4-176, the 
question of violation of this condition does not arise. 

Accurac y s The error bound as given by the relation 

2 

4- 196 is of the order 0( At Ax ), ana at times a i,7eaker error 

2 

estimate of the order 0( At + Ax ) is also reported [115]. 
Smaller values of Ax and At will without doubt give more 
accurate result biit the stability condition given by equation 

5- 5 does not permit infinitely small values of Ax. Thus the 
raaximuffi number of axial grids given by the relation 5-5 fixes 
the size of Ax, Therefore, qualitatively speaking, there 
does not seem to be much scope for the variation of Ax to 
Increase the accuracy. 

(iii) Eh.vsical and Boundary Conditions ; 1 typical 
cylindrical sample used in the experiment has a length of 4 ems, 
an inner diameter of 0.4 cm, and an outer diameter of 1,0 cr 
Prom physical considerations, ten axial grids over half the 
length of cylinder (i.e, 2,0 cms.) should give results that 
are fairly satisfactory, but. the maximum limit for the number 
of axial grids being six, a compromise is called for the use 
of axial grids less than six, say five (i.e. Ax = 0.2), On 
the other ha.nd, an attempt to reduce the number of grids below 
five will no doubt decrea,so the iiiaoh.in.Q tinie of coDip'u.ij^'tion. 



150 


but tile results are not likely to yield accurate wall profiles, 
Hie given boundary condition demands that the concentration 


at the pore mouth remains throughout the course of reaction 
In most cases, the reaction being diffusion controlled, the 
fall of pore mouth concentrations across the axial distance 
is very rapid. This would call for the use of a smaller size 
of axial grids and perhaps a different finite difference 
scheme near the pore mouth, thus leading to an increase in 
the computation time wliich is already on the higher side. The 
size of the grids, therefore, cannot be altered on the basis 
of the boundary conditions. 

Taking into consideration all the atove factors, A x = 0.2 
appears to be a judicious choice for the axial grid. How, 
using the Von Neumann condition 5-3, the maxiaium value of At 
can he chosen on the basis of b=0.5. The lov/er limit for At 
is decided by the magnitude of the round off errors in the 
computations, (see equation 4-197). It can be proved [119] 
tha^t for b=l/6, the error is of the order Ax^*^', and for all 
practical purposes, the corresponding value of At gives results 
vvhich are very close to the anal 3 :'tical solution if it exists. 
Thus the proper range of variation of At for the purpose of 
numerical trials appears to correspond to the variation of b 
between l/6 and 2/5. In the ultimate analysis, the most suited 
values of ax, At, and b are decided by a trial and error 
procedure so as to obtain maximum accuracy and minimum time of 
computation. The relative merits and demerits are illustrated 



151 


lay means of numerical data given in Table 5-2. 

It is clear from Tabic 5-2 that the difference between 
the maximum (0.238464 eras) and minimum (0.216568 eras) values 
of pore mouth radius calculated indeijendently by means of 
various finite difference schemes after half an hour of 
reaction at 1373°K is 9.18?^, which is fairly reasonable consi- 
dering the range of variation of Ax between 0.1 and 0.2 and 
the radius of convergence b betv/een l/6 and l/2. A comparison 
of the results of various difference schemes for the same mesh 
sizes ( = 0,2 and b = l/6) shows that the Modified . Rechardson 

Difference scheme gives the maximum value of 0.23 1001 eras, and 
the modified bacloward difference scheme gives the minimum value 
of 0.229986 eras. The variation is therefore 0 , 866 ^, Similarly, 
for Ax = 0,2 and b = 0.4, this variation is 1,43^. Presuming 
that the value b - l/6 yields results very close to the ana- 
lytical solution, the round-off and discretization errors 
involved in a numerical solution with b=0,4 may be estimated 
to be of the order of 5^u A similar trend is noticed in the 
results of the calculations for the reaction at 875°K. 

As regards the machine time [lEIT-7044]for computation, 
it can be easily seen from Table 5-2 that the Explicit Forward 
Difference schemo with Ax = 0,2 and b = 0,4 takes the minimum 
time of 835 seconds. This also includes the compilation time 
(of the order of few seconds) which however is negligible. 

It may be noted that convergence for the nonlinear simultaneous 
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R} =0.2 cm V 
Wo=4.25 gms 
At =0.06202 sec 


:; , 153 ^ 

k = 0-1545 cm /sec D = 1-032 cm^/sec 

Re=0-5cm H = 2.0cm5 

^Ao=0-292x10gmoles/cc Ax = 0-2 
Convergence factor =0.4 


t HOURS 

3.0 

28-0 

3-0 

(WqW^) g m 

0.058 

0-695 

0-023(exptl.) 



PROFILE DATA FOR 873'’ K 


FIGURE 5-6 




k = 501cm/sec 
Re=0 02 cm_^ 

C^q= 0-2 275x10 g moles/cm 

Convergence factor=0-4 


T =1123 K 
Rj = 0;C1cm^ 
VV^=Q-0162 gm 

At=0Q4211sec 


D=1-52cm/se 

B=2-0cms 

AX=0-2 


WALL PROFILES 


CONCENTRATION PROFILES 


PROFILE DATA FOR 1123 K 
figure 5-7 , 


t HOURS 

0-5 

10 

2-0 

3-23 

(W5W^)gnn 

4 2xicF 

97x1^ 

2-2xTd^ 

5-8xT(^ 
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T =1373 k=537cm/5ec. D = 2 04cnr?/sec. 

Ri = 0-2 cm. Re=0.5cm. H= 2 0cms. 

Wq= 4-25 gms. C^Q=0-186x10gmoles/cG Ax=0’2 

At=0-03156sec. Radius of covergence -0-4 



PROFILE DATA FOR 1373 K 
FIGURE 5-8 : 




T=1573 K 
Ri=0-01cm, 
Wo=0.C162 gm. 
Atr0-0257sec. 


V - : \ 156 

k = 68-3 cm/sec. D =2-49 cm/sec . 

Re=C-02 cm. H = 2-Ocms. 

^ " ■ -5 . 

CAo=0.1623x10g.moles/cc AX40.2 

Convergence factor = G-4 


;t HOURS ; 

0.50 

TO 

1-50 

1.8 

(Wo-Vv.(-) gm 

0-79x10^ 

1.835x10"^ 

3 Sxio'^ 

5.4x10^ 




CONCENTRATION PROFILES 


PROFILE DATA FOR 1573 K 








1b/ 


T= 2073 k = 229,cnn/sec. 

Rj =0r01 Re=0'02cnn H=2-0cms. 

V\^ = 0.0162 g ms. CAo=0.1232xic?gmoles/:c AX = 0-2 

At =0-0l698 sec Convergence factor=0-4 



Concentration was almost zero throughout 
the length. At the pore mouth C = *0001 


2073"'' k 


FIGURE 5-^10 
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T- hypothetical 
Rj = 0.2Gm 
Wo = 4.25 grms 
At -G.032 sec. 


GAo=1’595x10gnnoles/c.c. AX=0-2 
Radius of convergence-0-4 


WALL PROFILES 


CONCENTRATIQN PROFILES 


FIGURE 5-11 
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algebraic equations in Rose’s sclieine was not persued because 
of tne prohibitively long time of computation and probably 
no added advantage in accuracy. Taking into consideration 
the above factors computations were done on IBM 7044 and 
IBM-360 using the Explicit Eor'ii'ard Difference gcheme for all 
temperatures of reaction with Ax = 0.2 and b = 0.4. The ' 
concentration and wall profile data thus obtained are graphi- 
cally represented on pages 153 - 158. 

It is clear from the plots on pages 153-1 58 that the 
numerical data for axial concentration distributions and pore 
wall profiles are too diverse for quantitative comprehension. 
The length of computation time is the most crucial factor 
which puts serious limitations on the accuracy and vastness 
of the data, or for that matter, on the general scope for 
qualitative or quantitative betterment of the data. With a 
view? to rediice the computation time, several alternatives in 
Eortran language w?ere tried. Accuracy, being very important 
in the initial stages of computation, the variable C (of the 
order of 1) was replaced by a nav variable f = 1-C (of the 
order of zero) so that the same 16 significant decimal places 
retained by the computer gave greater accuracy. Besides 
attempting various difference schemeg the distance and time 
grids were also increased to their maximum permissible values. 

In this connection, it may be noted that the five percent error 
introduced by the increased time increment At is alw/ays on 
the negative side i.e. the poremouth radius at hi^ier At is 
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always less tiian that at lower At, Although this error 
affects the quantitative results, the nature of the quali- 
tative conclusions remains unchanged. Thus, unto the time 
an analj/tical solution or a quickly converging finite 
difference scheme for the system of quasilinear parabolic 
differential equations describing the unsteady state chemical 
reaction in a pore are made available, the conclusions drawn 
from single pore stijdies may at best be considered, qualitative 
comprehension. 

Before going into the details of conclusions, some 

observations may be made regarding the experimental anomalies 

and computational procedures. Graphs on. pages 81 -.and 82 indicate 
■bhat higher temperatures fail to give 4n.cr.ea^ed rate constants. 


In fact, if experim.ental accuracy is taken for granted, k at 
870° C actually is less than that at 746° C, This may apparently 
look absurd, but if Knudsen diffusivity is operative, the 
results become well explainable and the levelling off of 
Vireight loss curves is quite a common phenomenon under the 
•diffusion controlled regime. Secondly, on page 155> the 
profile data indicates that the experimentall^r observed pore 


mouth radius is less than that predicted by computations, 

A similar trend is noticed at other temperatures, too. This 
general deviation may be due to a number of reasons. 


Por beginning calculations in the Modified Richardson 
Difference Scheme, independent values were chosen from ■3J§:BUlts oj 
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explicit fomA'ard difference scheme and Milton lees Extra- 
polated Scheme. Eyen in other schemes, to doubly ensure 
stcLi^ibi, values obtained from one difference scheme were 
fed as the starting data for another scheme, and the error 
trend was observed to be unchanged. In most of the cases, ' 
computations were stopped when the poremouth radius became 
equal to the outer radius because by then the quasisteady 
state gets well established. In cases v/hore calculations 
are continued beyond this stage, there is no solid material 
left for reaction at tne pore mouth. Therefore equation 4-122 
IS replaced by the corresponding diffusion equation. 

f(n,m+l) -f(n,m) = Tf(n+i ,ni)-2f (n,m)+f(n-1,m)| ...5-6 

and R(n,m+1)= R(ntm) ...5-7 

Thus forward march is continued by using a suitable combination 
of equations 5-6, 5-7 and the respective difference forms 
for 4-122 and 4-123, Mesh sizes in eqa.5r'6 are Jcepl? unaltered 


in order to retain the same magnitude of error. This is one 


more requirement which dictates 
Difference schemes mentioned so 



the choice of mesh sizes, 
far give the radius only at 
grid points 1,2, 3, 4 etc. (see 
figure 5-2), There is no 
provision to calculate the 
rate at which point A travels 
towards point D at grid number 2. 
It is difficult to calculate this 
without increasing the number 


EIGURE 5-2 
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of grids ax the pore mouth and hence i/i-ithont affecting the 
stability and accuracy in the course of calculations. To 
avoid this, the average rate at which A travels towards D 
is taken as the rate at which C travels towards B, Thus the 
rate at which the interface travels is determined. 

As has already been pointed out in the introductory 
remarks, the use of the parameter ’temperature- in improving 
the yield of a heterogeneous gas solid reaction-‘is generally 
overstressed. To bring to surface the fact that temperature 
can not only be an unimportant factor, but even a negative 
factor in improving the yield of a chemical reaction of a 
porous solid, a combination of physically realizable k, D and 
^Ao their variation with temperature is presented in 
Table 5-3. 


TABLE 5-3 


Temperature 

■ 

cm 

s ec 

2 

sec 

Ao cc 

Wt. loss of a stan- 
dard cylinder after 
4- hour reaction. 

b 

5 

1.50 

0.23 X 107^ 

13 mg 

m 

xp 

'50 

1.75 

0.19 x 107^ 

22 mg 

m 

"3 

■ > i 

- j 

100 

2.00 

- 

0.14 X 10“^ 

18 mg 

Let 

Ic = 

e-E/HT 


...5-8 


h “ T“ ...5-9 

Cao=^ T“P ...5-10 

where k^ is the pre-exponential factor and E is the activation 
energy, E the ideal gas constant, and a and p are positive 







163 


constants, Now, if the reaction is diffusion controlled, and 
further, if the quasisteady state has already been attained, 
the ra.te of chemical reaction (H) is essentially governed by 


Ci'v 


the term or ~~ , where the factor ^ (pre- 

dominently dependent on pore geometry) remains constant with 


temperature, therefore R' 


DC , 

AO 

' r' T 
a-p. 


a „ “p 






Under quasisteady state conditions, the rate of diffusion 
controlled chemical reaction in a pore 

(a) increases marginally with temperature if a> P 

(b) does not change with temperature if a = p 

(c) decreases marginally with temperature if a p . 

For most of the systems at moderate pressures, p 1 and 0,4<a 
<2,2, If En.y,<isQi3;: diffusion is operative, =0.5 and it can 

be seen that the reaction rate actually decreases v/ith increasing 
temperature even though k increases. Maahs [lOO] identifies 
a similar phenomenon by referring to the deviation of k from 
normal Arrhenius plot for the oxidation of pyrolytic graphite 
at high temperatures, and keeps the discussion inconclusive 
by making a casual reference to the possible shortcomings of 
Scala’s . theory [101, 102] of diffusion controlled chemical 
react long. An accurate estimate of the index a for the variation 


of effective diffusivity with temperature may give a clue for 
the deviation of k, especially when the material is almost 
nonporous or the average pore radius is sniall. 
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Ilormally, for a diffusion controlled reaction the 
chemical reaction rate constant is sufficiently large and 
the Eiiolal concentration of the fluid reactant on the pore 
surface is nearly zero. Therefore the rate should be governed 
essentially by the mass transfer rate at the pore mouth. 
Experimental veigiit loss curves for most of the porous solids 
indicate a varying rate of reaction in the initial stage 
and a constant rate in the final stage. This phenomenon is 
easily explained on the basis of a "shrinking pore model", 

Eigure 5-5 illustrates a typical diffusion controlled fluid- 
solid reaction. The natural pore system indicated in 5-51 
is modelled by a parallel pore system shown in 5-5B, Bow, 
during the course of reaction the solid undergoes a surface ti: 
transformation from5-5B to 5-5C, and this results in a change 
in the po rewall profile and the concentration gradient at the 
pore mouth. This naturally causes a variation in the effective 
mass transfer coefficient at the interface. As the poremouth 
radii of two adjacent pores increase to their maximum, as shown 
in Eigure 5-3C, the quasisteady state sets in and the- interface 
recedes at a uniform velocity giving rise to a linear weight loss 
curve. The graph on page 81 clearly indicates a varying ^ 
followed by an almost constant dW/dt. This behaviour may be 
more pronounced when the pore-mouth concentration gradient 
is held constant by allo?i/ing the flow right across the pore 
mouth. In an actual porous solid, even though the pore geometry 
is very different from that modelled in Eig. 5-3B, weight loss 
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FIC-URB 5-3 

curves follow the above ursteadystate - quasi-steady state 
pattern with surprising a.ccuracy. A large number of pores 
?/]'iich statistically average out the deviations arising out of 
the modelled pore structure, may well explain this resemblance 
of behaviour between the actual and the modelled pore structures, 
Hovvever, this qualitative resembelencc cannot and should not 
induce an Estimation of the difficulties involved in the 

quantitative matching of theory and practice. At this stage, 
it may be worth noting that the so called diffusion controlled 
regime actnaally refers to the quasisteady state portion of the 
general diffusion controlled reaction. The unsteady state 
portion is generally omitted from the discussions or is mistaken 
for a chemical reaction controlled regime. This may be one 
of the reasons for the existing contraversies regarding transition 
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zone between the chemical reaction controlled, regime and the 
diffusion controlled regime. 

In the light of the above discussions, it may be 
interesting to explore the common belief that high density 
and low temperatures are desirable parameters for the use of 
graphite as a material of construction in rocket nozzles. It 
is clear that slightly porous material with fine pore distri- 
bution corrodes at a rate slower than high density material 
(nonporous solid). It has also been seen that when Kniidsen 
diffusion is operative, reaction rate decreases at high temp- 
eratures., Thus, at least from the point of view of chemical 
corrosion, hi^ly nonporous material and low temperature of 
operation are not alv^ays necessarily desirable. 

Coming to the 'jse of single pore data to calculate the 
effectiveness factor, it may be noticed that the poremouth 
concentration cannot always be used as the basis for the 
definition, hecauae the pore mouth itself changes its position 
during the course of reaction. The mouth concentration is 
thus affected by the external diffusion process. Therefore, 
to avoid confusion and unnecessary mixing up of initial pore 
dimensions v/ith the effectiveness factor, the effectiveness 
factor may be defined only for the period during Yihlch the 
inner radius increases and becomes equal to the outer radius 
at the poremouth, and only for that axial distance over which 
the inner radius actually varies. Variation of internal surface 
area with time becomes an added factor complicating the rate 
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expression over the entire surface. In view of these diffi- 
culties and shortage, of computation time, only an approximate 
picture of variation of effectiveness factor with time is 
indicated in figure 5-4. It is to note that these plots 
based on single pore. data have a remarkable similarity to those 
defined for entire porous solids [3, 6 ]. 



flG-UZB 5-4 . FIGURE 5-5 


Similarly in Figure 5-5j the variation of concentration 
and the rate of weight loss are qualitatively indicated. Ini- 
tially the rate of weight loss decreases and then it gradually 
increases as the quasisteady stale approaches. In a typically 
diffusion controlled reaction the percentage conversion during 
the first half hour is nearly half of that during the subsequent 
half hour. This fact will have an important effect on the 
adjustment of the residence time of a pellet in the reactor. 

It may be profitable to maintain the reactor at less costly 
conditions in the initial stages and more costly conditions 
during later stages to attain maximum conversion at minimum cost 
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More quantitative data and detailed investigations are called 
for to impliment such schemes in practice, 

from the concentration profiles plotted on pages 155jl56, 
157s it is clear that in many cases a major portion of the 
interior of the pore has a fluid concentration which is almost 
zero. Thus an effective utilization of huge pore surface area 
would he possible only when high diffusivities are attained. 

An obvious method would be to introduce pressure diffusion. 

While empirical data regarding the effect of pressure on the 
rate of production in blast furnaces, coke ovens, and limekilns 
etc are available, satisfactory mathematical models have not 
heen proposed. The reasons for the same may be the occurrence 
of the bulk, flow term in the material balance equations, the 
varying order of reaction and the dependence of ciffusivity on 
concentration at high pressures. 

Since it is the first time an attempt is made to 
develop a '’shrinking pore model" v^ith an objective of defining 
precisely the quantitative effects of the physical factors on 
the chemical reactivity of porous solids, and of deciding the 
optimal conditions for carrying out heterogeneous reactions, many 
difficulties are likely to be experienced in modifying or 
sophisticating the present model on the following lines: 

(a) Introducing radial diffusivity term in the model, 

(b) Defining the first order rate constant on the 
volumetric basis, and introducing it as a 
function of inner radius. 
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(c) Considering the axial diffusion in gas phase 
and radial diffusion in solid phase, to set 
up_^two different partial differential equations 
and solving them simultaneouslyi 

(d) Express k, D, and as functions of temperature 
in the present model^and try to obtain an 
expression for the optimum temperature for the 
heterogeneous gas solid reaction on the basis of 
"shrinking pore model," 

(e) Use the present numerical data to get an analytical 
expression for the pore radius in terms of time 
and axial distance, and introduce this R(x,t) in 
the differential equation 4‘-l22, and solve it 
again to get a better approximation R,,(x,t), and 
continue the iteration until accurate an^ y tic al 
expressions for R(x,t) and C(x,t) are available. 

(f) Express k as a fi'inctionof time varying pore 
struct'are and incorporate it into the main 
equation 4“122. 

(g) Use the model for different orders of reactions 
or for series and parallel reactions under non- 
iso thermal conditions, 

(h) Y/hen there is a change in the number of moles 
of the gaseous reactant during the course of 
reaction, introduce a term for bulk flov/ due to 
pressure or express I) as a function of concen- 
tration and solve the resulting nonlinearity, 

(i) Introduce a separate terra for pressure diffusion, 

A large number of such modifications are possible, 

but the key factor for tho success of any such modification 
is the availability of a suitable difference scheme which 
yields satisfactory results by using a reasonably short time 
of computation, A number of new difference forms [112] are 
being reported and it is hoped that the v;hole range ot design 
data for a heterogeneous system - v/hich till now are inadequate 
to cope up Yi/ith the complex problems encountered even in rela- 
tively simple processes - world soon be available. 
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APPENDIX - A 

One of the common methods of linearizing the equation 4-3 
for the determination of k from the experimental weight loss 
curves is to make an approximation of. the following type 

2 -“I ■ - 

e^ {T-erf (x)) = [ ax + b] ...A-l 

where a and h are nuaierical constants and x = k - y /t. 

The approximate relation A-1 can be made accurate at two values 

of X by choosing suitable numerical values of a and b. Assuming 

maximum and minimum limits for x as 2.0 and 0.05, the equation A-' 

may be rewritten as 

2 , - ^ ■ - -1 
e^ {1-erf (x)} = [ 1 .772x +’ 0.736 ] ...A-2 

Thus, equation 4-30 becomes 


1 ^ 

S dt 


= k C 


[ 1.772y /t + 0.736 ] 


-1 


...A-5 


0.736 ,dW^ -1.772y / dW ^ 
S ^dt^ - S H ^ 

dW r-r* dW 


k C 


A. 


.. .A-4 


\jL Vv M V i 

Thus a plot of ^ versus / x should result in a straight 

line, the intercept of which gives k C. and hence the value of 

% 

k. It may be noted that the maximum error obtained by differen- 

2 V --I 

tiating the function e^ {1 - erf(x)} - [l.772x+ 0.736 3 v/ith 
respect to x is of the order of 45^. Thus simplifications of 
the type A-2 become useful in the estimation of the first order 
rate constant. 



